code: def TEME_to_ITRF(jd_ut1, rTEME, vTEME, xp=0.0, yp=0.0):
      
      theta, theta_dot = theta_GMST1982(jd_ut1)
      zero = theta_dot * 0.0
      angular_velocity = array([zero, zero, -theta_dot])
      R = rot_z(-theta)
  
      if len(rTEME.shape) == 1:
          rPEF = (R).dot(rTEME)
          vPEF = (R).dot(vTEME) + cross(angular_velocity, rPEF)
      else:
          rPEF = einsum('ij...,j...->i...', R, rTEME)
          vPEF = einsum('ij...,j...->i...', R, vTEME) + cross(
              angular_velocity, rPEF, 0, 0).T
  
      if xp == 0.0 and yp == 0.0:
          rITRF = rPEF
          vITRF = vPEF
      else:
          W = (rot_x(yp)).dot(rot_y(xp))
          rITRF = (W).dot(rPEF)
          vITRF = (W).dot(vPEF)
      return rITRF, vITRF
docs: Convert TEME position and velocity into standard ITRS coordinates.
  
      This converts a position and velocity vector in the idiosyncratic
      True Equator Mean Equinox (TEME) frame of reference used by the SGP4
      theory into vectors into the more standard ITRS frame of reference.
      The velocity should be provided in units per day, not per second.
  
      From AIAA 2006-6753 Appendix C.