Need help for GMT x2sys adjustment

Hi, Prof.Paul
I guess this is a question for you (creator of the X2SYS).

I am attempting to do the crossover adjustment of the satellite altimetry sea surface height data using the x2sys. I learned the processing through your example codes of paper Wessel, P. (2010), Tools for analyzing intersecting tracks: the x2sys package, Computers & Geosciences, 36, 348-354.

The determination of the crossover locations and the SSH differences using x2sys_cross works well.

However, the adjustment of the crossover differences using x2sys_solve seems not right. And finally, after adding the correction using x2sys_datalist, I see a wrong figure.

I think the x2sys should work for satellite altimetry SSH adjustment since the aim of all data adjustments is similar.

Here is all the data and GMT (748.9 KB)

Looking forward to hear you.


Just posting the code here from the zip file:

#!/usr/bin/env bash

export X2SYS_HOME=.
rm -rf CROSSHY2

(cd hy2; ls *.dat) > hy2.lis

gmt x2sys_init CROSSHY2 -Dgeoz -Edat  -G -F 
echo "hy2" >> CROSSHY2/CROSSHY2_paths.txt
gmt x2sys_cross -TCROSSHY2  =hy2.lis -Qe -R140/200/0/50 -Ia -V > outs.txt
gmt x2sys_list -TCROSSHY2 outs.txt -Cz -Fndc -Vl > COE_use.txt
gmt x2sys_solve -TCROSSHY2 COE_use.txt -Cz -Ed -Vl> corrections.txt
gmt x2sys_datalist -TCROSSHY2 -Lcorrections.txt =hy2.lis -Flon,lat,z >tmp.dd
# plot
gmt blockmean tmp.dd $R -I10m | gmt surface $R -I10m -T0.25
gmt grd2cpt  -Crainbow >t.cpt
gmt grdimage $R -JM8c -K -Ct.cpt >
gmt psbasemap -R -J -P -K -Ba -O >>
gmt pscoast -R -J -W1/0.25p -Dl -Glightyellow -P -A1000 -O -K >>
gmt psscale -Ct.cpt -D2.0i/-0.5i+w8c/0.15i+h+jTC+e -O -Bxaf -By+l"m" >>

#Obtain adjustments and grid the corrected and adjusted data
gmt x2sys_report -TCROSSHY2 -Cz outs.txt -Lcorrections.txt -A > report.txt

gmt psconvert -A -P -Tg

I’m not that great an x2sys expert like Paul, but I do work a bit of work on satellite altimetry (cryosphere rather than ocean surface) and have worked on the PyGMT x2sys_cross wrapper (which you’re welcome to try :grin:). Without having looked too closely into your code/data, here’s some things I would try:

  • Use -D to turn off the cylindrical polar conversion, it looks as though these are latitude based errors?
  • Different interpolation (-I) settings. I noticed you’re using Akima spline (-Ia), maybe try linear (-Il) or cubic (-Ic) to see if it makes a difference
  • Use more data points -W than the default -W3 and see if it makes a difference.

That said, maybe it’s just your colorbar being binned to equal-ish intervals. Try using continuous with -Z and see if it helps.

Thanks :wink:

At present, x2sys_cross looks fine to calculate the ssh difference at crossover points. So -D is probably not the reason. Also, -I and -W will probably affect little on the result.

We do know there may be an issue with x2sys_solve. Have you tried the simplest correction first, i.e., -Ec ? I suspect the problem has to do with the linear drift rate calculations or distances.

I just tried the -Ec. It is also wrong.

This is the original SSH:

We can see the sea surface height is changed with a big scale after adjustment.

I played with just doing the constant shifts. As you know, getting crazy offsets ranging from -9000 to 6000. All offsets sum to 0, which is the constraint we give it. Could you please help me debug this by identifying 3 descending and 3 ascending tracks that are next to each other, so that we should find 9 crossovers in the middle of your area. I.e., find which tracks those are. It is too hard to look at a 156x156 matrix… Just post that updated hy2.lis with those 6 tracks and I can have another look.

Oops, one more update. When I only solve for constant corrections it works fine. I do this:

gmt x2sys_list -TCROSSHY2 outs.txt -Cz -Fnc -Vl > COE_use.txt
gmt x2sys_solve -TCROSSHY2 COE_use.txt -Cz -Ec -Vl> corrections.txt

and the rest of the script is the same. Because the corrections are pretty small ( -0.322769/0.654737, which they should be) the grid2.png plot is just a little bit different than grid.png. So when you said it did not work, perhaps you did not use the options I did? You would need to change both list and solve.

The linear drift solution is still crazy so still need those 6 tracks identified.

Hi,Prof. Paul
Thanks very much.

Here is the 6 (68.2 KB)

Yes. It works fine without the distance involved.

Not sure what is going on. 13x13 matrix problem for these 6 tracks is very unstable. The simple test example in GMT also uses distances but yields reasonable results so I don’t think the code is wrong. I can tell in debug that the individual adjusted crossovers (after both tracks are adjusted) shrink, hence the mean/std reported improves, but the values of the adjustments are huge. The just happen to cancel out mostly when computing crossings. For just constant shifts things are stable. I think I need to add another condition on the solution. Sorry, this wont be solvable quickly. I suggest you stick with just the constant adustment for now.

Thanks again.
It is OK. This is not urgent for me. I just find the x2sys is an interesting and powerful tool and could be used in the satellite altimetry processing.
We also have other tools to do crossings and adjustment. But slow than GMT.

1 Like

I seem to remember that when solving for COEs in magnetic anomalies the results were crazy when cross-overs happened very far from the beginning of the lines.

The old XSYSTEM in GMT 4 solved this problem by iterations and it was stable. It only knew of offset plus slope * time_since_start_of_track. x2sys generalized all this and formed a single set of normal equations. However, I think I understimated (meaning, I was a bit clueless) about how to ensure the constants don’t fly all over. I guess I was fooled by the initial test cases that worked fine. Will check with Leo on that linear system for advice.

My gut feeling is: if the angular coefficient is wild but the cross-over conditions are met, then this is a sign of instability and/or non-uniqueness in the adjustment. To fix this, you’d need to introduce constraints:

  1. New observational conditions for the parameters. For example, the predicted data should have a known value in some places or conditions like closed loops.
  2. Regularization. In this case, some damping would go a long way. But then you have to be more careful with the scaling of the objective function to avoid crazy damping parameter values.

It would be much better to try the first solution before going for regularization since that gets more subjective. I’m not very familiar with the methodology in x2sys so I’ll have a look at Paul’s paper to see if anything comes to mind. I feel like it should be possible to add extra conditions to better constrain the problem.

Hi yangleir, may i know what is another tools to perform crossover adjustment?

The gravsoft could do adjustments of altimeter data. But it is not open source.

Hello I’m new here,I need help in doing crossover adjustment in satellite altimetry sea surface from my undergraduate thesis. I have already downloaded the data sets from RADS. what software do I need to perform this crossover adjustment?

Suggest you start reading this paper first, then begin with x2sys_init.

1 Like