I spent a few more hours looking at this but do not find anything wrong in our implementation.
I added the following statements in the test code:
Rotation r0 = frameTOE.getTransformTo(FramesFactory.getEME2000(), date).getRotation();
Rotation r1 = frameTOD.getTransformTo(FramesFactory.getEME2000(), epoch).getRotation();
Rotation r2 = frameTOD.getTransformTo(FramesFactory.getEME2000(), date).getRotation();
System.out.println(date);
System.out.format(Locale.US, "r0 TOE -> EME2000: %18.15f %18.15f %18.15f %18.15f%n",
r0.getQ0(), r0.getQ1(), r0.getQ2(), r0.getQ3());
System.out.format(Locale.US, "r1 TOD(epoch) -> EME2000: %18.15f %18.15f %18.15f %18.15f%n",
r1.getQ0(), r1.getQ1(), r1.getQ2(), r1.getQ3());
System.out.format(Locale.US, "r2 TOD(date) -> EME2000: %18.15f %18.15f %18.15f %18.15f%n",
r2.getQ0(), r2.getQ1(), r2.getQ2(), r2.getQ3());
System.out.format(Locale.US, "P TOE: %14.4f %14.4f %14.4f%n",
pvTOE.getPosition().getX(), pvTOE.getPosition().getY(), pvTOE.getPosition().getZ());
System.out.format(Locale.US, "P TOD = TOE.getTransform(TOD) (p TOE): %14.4f %14.4f %14.4f%n",
pvTOD.getPosition().getX(), pvTOD.getPosition().getY(), pvTOD.getPosition().getZ());
Vector3D pEME2000 = r0.applyTo(pvTOE.getPosition());
System.out.format(Locale.US, "P EME2000 = r0(P TOE): %14.4f %14.4f %14.4f%n",
pEME2000.getX(), pEME2000.getY(), pEME2000.getZ());
Vector3D pTODBis = r2.applyInverseTo(pEME2000);
System.out.format(Locale.US, "P TOD bis = r2^-1(P EME2000): %14.4f %14.4f %14.4f%n",
pTODBis.getX(), pTODBis.getY(), pTODBis.getZ());
These allow to track down to the most basic thing (3D rotations only), your test cases. Hare are the results I get, with IERS 2010 conventions and EOP from the latest finals20001.all file from IERS:
2008-11-15T08:30:00.000
r0 TOE -> EME2000: -0.999999386362698 -0.000014694334147 0.000441454418156 -0.001015960776932
r1 TOD(epoch) -> EME2000: -0.999999386362698 -0.000014694334147 0.000441454418156 -0.001015960776932
r2 TOD(date) -> EME2000: -0.999999386050600 -0.000014659565422 0.000441566711815 -0.001016219641631
P TOE: -34810228.1700 -23815213.8700 -131484.6301
P TOD = TOE.getTransform(TOD) (p TOE): -34810215.8114 -23815231.9002 -131490.7947
P EME2000 = r0(P TOE): -34858649.1421 -23744436.3818 -100030.1197
P TOD bis = r2^-1(P EME2000): -34810215.8114 -23815231.9002 -131490.7947
and
2013-07-26T19:25:00.000
r0 TOE -> EME2000: -0.999998579469317 0.000016809332261 0.000671638270993 -0.001545858605100
r1 TOD(epoch) -> EME2000: -0.999998579469317 0.000016809332261 0.000671638270993 -0.001545858605100
r2 TOD(date) -> EME2000: -0.999998579313809 0.000016780989368 0.000671675324374 -0.001545943407975
P TOE: 11060027.4100 -365735.6600 -675152.5700
P TOD = TOE.getTransform(TOD) (p TOE): 11060027.5220 -365733.7454 -675151.7721
P EME2000 = r0(P TOE): 11057926.9337 -399903.9909 -690020.7256
P TOD bis = r2^-1(P EME2000): 11060027.5220 -365733.7454 -675151.7721
So I extracted the pure rotations between TOE and EME2000 (r0), and between TOD and EME2000 at epoch (r1) and at date (r2). As expected, r0 and r1 are exactly the same, and r2 is slightly different (differential effect of precession and nutation during the day).
Then I converted the position, first using the complete TOE.getTransformTo(TOD)
and then simply applying rotation r0 to the TOE position in order to get the position in EM2000 and then applying the inverse of rotation r2 to the position in EME2000 to retrieve position in TOD. I still get the exact same results.
I would advise you to check with your reference program the quaternions with respect to EM2000 and the EME2000 position also.
Orekit frames have been extensively validated by different teams and are known to be extremely accurate (centimeter level or better at Earth radius). The tests with the rotations above make me think the getFrozenFrame()
underlying algorithm does what it is expected to do: it freezes the TOD/EME2000 rotation to a constant which corresponds to the freezing date, and then it applies it throughout timeline. I don’t see any errors here.