code: def trans_coeff(eq, x, y, z):
      
      tt = (eq - 2000.0) / 100.0
  
      zeta = 2306.2181 * tt + 0.30188 * tt * tt + 0.017998 * tt * tt * tt
      zetto = 2306.2181 * tt + 1.09468 * tt * tt + 0.018203 * tt * tt * tt
      theta = 2004.3109 * tt - 0.42665 * tt * tt - 0.041833 * tt * tt * tt
  
      zeta = math.radians(zeta) / 3600.0
      zetto = math.radians(zetto) / 3600.0
      theta = math.radians(theta) / 3600.0
  
      p11 = (math.cos(zeta) * math.cos(theta) * math.cos(zetto) -
             math.sin(zeta) * math.sin(zetto))
      p12 = (-math.sin(zeta) * math.cos(theta) * math.cos(zetto) -
             math.cos(zeta) * math.sin(zetto))
      p13 = -math.sin(theta) * math.cos(zetto)
      p21 = (math.cos(zeta) * math.cos(theta) * math.sin(zetto) +
             math.sin(zeta) * math.cos(zetto))
      p22 = (-math.sin(zeta) * math.cos(theta) * math.sin(zetto) +
             math.cos(zeta) * math.cos(zetto))
      p23 = -math.sin(theta) * math.sin(zetto)
      p31 = math.cos(zeta) * math.sin(theta)
      p32 = -math.sin(zeta) * math.sin(theta)
      p33 = math.cos(theta)
  
      return (p11, p12, p13, p21, p22, p23, p31, p32, p33)
docs: This function is provided by MOKA2 Development Team (1996.xx.xx)
      and used in SOSS system.