Question about the generation of the nonConservativeForces and neutralDensity products

[Levi-always-rejoice]

Dear GROOPS developers,

I am currently working on the retrieval of non-conservative accelerations and thermospheric density from precise orbit products of GRACE-FO, and I have a few questions regarding the operational products provided by TU Graz.
I noticed that the operational dataset includes both the nonConservativeForces product and the neutralDensity product. However, I could not find any documentation or publication describing how these products are generated.

My first question is:

1. How are the nonConservativeForces and neutralDensity products computed?

• Are the non-conservative accelerations obtained by applying Lagrange interpolation and numerical differentiation to the precise satellite positions or velocities, similar to the approach proposed by Calabia et al. (2015, https://doi.org/10.1016/j.ast.2015.05.013)?
• Or are they estimated directly during the reduced-dynamic POD process?
• Is there any documentation or publication describing the processing strategy in detail?

I also performed my own computation using the JPL 1 Hz reduced-dynamic orbit products.
Specifically, I

• computed the total acceleration from the precise satellite velocities using Lagrange interpolation and numerical differentiation as introduced in Calabia et al. (2015, https://doi.org/10.1016/j.ast.2015.05.013),
• computed the conservative accelerations with GROOPS (including gravity field, solid Earth tides, astronomical tides, Doodson harmonic tides, pole tide, ocean pole tide, and relativistic effects),
• and obtained the non-conservative acceleration by subtracting the conservative acceleration from the total acceleration.

When comparing my non-conservative accelerations derived from the JPL 1 Hz reduced-dynamic POD product for GRACE-FO-1 on January 1, 2023, with the corresponding TU Graz nonConservativeForces product, I found that the overall temporal variations in the along-track component agree reasonably well, but there is a noticeable bias between the two. (see figure below).

Interestingly, after converting the non-conservative acceleration into thermospheric density, my retrieved density agrees very well with the TU Graz neutralDensity product in terms of the overall temporal variations. However, my retrieved density contains noticeably stronger short-term fluctuations than the TU Graz product. At this stage, I have not figured out whether these fluctuations originate from the POD products, the numerical differentiation, or some differences in the subsequent processing. (see figure below)

Therefore, I have two additional questions:

2. In the TU Graz nonConservativeForces product, the Y-component is approximately two orders of magnitude larger than the X- and Z-components (see figure below). In contrast, all three components in my results are of the order of 10……-7 m/s……2. Is this expected?

3. What could be the possible reason for the bias in the non-conservative acceleration while the retrieved densities still agree well?

• Could this be related to different POD strategies (e.g., reduced-dynamic vs. kinematic)?
• Or are there additional empirical accelerations, bias corrections, or other processing steps included in the TU Graz products?

Any clarification or references would be greatly appreciated.
Thank you very much for your time and for providing these excellent datasets and the GROOPS software.

Best regards,
Wei Li

[ajs8041]

Dear Wei,

Concerning your main question, there is currently no paper that describes our GRACE-FO density product in detail. However, I am in the process of writing one and I can obviously provide some input.

Most importantly, the GRACE-FO neutral density product on our ftp is based on accelerometer measurements and not on POD observations. However, POD observations are used to calibrate these accelerometer measurements (i.e., bias + scale are estimated in a least squares adjustment such that the ACT1B product fits to a POD solution). Thus, your purely POD derived densities should be a good fit and that seems to be the case.

Now to your pictures:
Concerning the first and third picture, according to my understanding, the nonConservativeForces product on our ftp for that time is actually uncalibrated accelerometer data. The biases for x and y look like this:


Since these values are considered for the calibration of the ACT1B product and the densities match well in your second illustration.

EDIT:
Why does this matter? If you check the nonConservativeForces product for a more recent date, .e.g., 2026-01-01, you will see that the large accelerations in y-direction are gone. This is because we usually simulate the ACC data (radiation pressure + drag as provided by DTM2020) at this processing step for availability reasons. I can not pinpoint the reason why this is not the case for 2023-01-01 and maybe other days. There will still be a bias (primarily of the DTM2020 model vs reality/POD), but this one is going to be relatively small compared to the one from 2023-01-01 you showed in your illustrations. So concerning the nonConservativeForces data product: at this period, we apparently used the ACT1B product (which contains a bias) but usually we use simulated data (which contains a smaller bias). For the when and why I am not familiar enough with our POD processing.

Maybe a final hint, since we also use kinematic orbits for processing densities. The GNSS observations are inherently noisy and in some sense they are much noisier than the accelerometer measurements. This should be reflected in the density products, POD derived densities need some smoothing strategy as frequencies below minutes tend to be noise and not much else. This is how our solutions for GRACE-FO compare:

Feel free to keep me updated on your progress, you are definitely on a good path.

Best regards,
Andreas

[Levi-always-rejoice]

Dear Andreas,

Thank you very much for your detailed and thoughtful reply. It has clarified many of the questions that I had, especially regarding the relationship between the nonConservativeForces, ACT1B, and neutralDensity products. I really appreciate you taking the time to explain these details.

I have one additional question out of curiosity.
In your last figure comparing the ACC-based and POD-based densities, the POD-derived density appears to be very smooth. However, the GRACE-FO-1 kinematic orbit product available on the TU Graz FTP is sampled at 10-second intervals and, as far as I can tell, contains only satellite positions.
I am therefore wondering how such a smooth and high-quality POD-based density can be retrieved from those orbit products. Is there an intermediate processing step (for example, a higher-rate orbit representation, orbit interpolation, smoothing, or another strategy) before deriving the density?

Of course, if this involves methods that will be described in your upcoming paper or information that is not yet ready to be shared publicly, please do not feel any pressure to answer this question. I completely understand and will be happy to read the paper once it becomes available.
Thank you again for your help. Your explanations have been extremely valuable and have given me much more confidence that I am following the right direction.

Best regards,
Wei Li

[ajs8041]

Dear Wei,

our processing of POD based densities does in fact use the 10s orbit data that you mention. I guess it does not matter too much at what point and how exactly you smooth the data, but I will outline how we do it anyway.

We use a somewhat different approach to the one from A. Calabia you mention in your original question. I will thus describe our smoothing strategy in a maybe heuristic way which hopefully is still helpful to you. To my understanding, you calculate an acceleration for every single data point without, on this level, linking the epochs together in any way. This means, any outlier is going to end up as an acceleration. On the contrary, we model the orbit without drag and then fit it to kinematic observations. And this is exactly where the smoothing happens, we do not fit single data points independently but rather link them together (the choice of the window length is very empirical, at the moment we use 15 minutes) and downweight outliers iteratively. This way, the accelerations are smooth.

I am not sure how many details I am going to include in my paper, since it includes several parameter choices that I would have to substantiate.. Also, the choice of how and when you smooth the data is probably not critical. I simply believe that the outliers in the GNSS/orbit data should be dealt with somehow as they are likely not physically meaningful.

Best regards,
Andreas