ldeo-glaciology / xapres

package for processing ApRES data using xarray
MIT License
3 stars 2 forks source link

Adding to_zarr bug fixes + preliminary range_diff calculations #14

Closed glugeorge closed 1 year ago

glugeorge commented 1 year ago

Since I cannot create issues on my forks, I will use this to track my changes and give updates. The first things that should be merged are the edits to the zarr notebooks, since that amends the issue of encoding complex numbers as floats in addition to the wrong date encoding for A103. I also have preliminary iterations on a range_diff function. At first, it could only calculate between two bursts, but I've spent the past couple of days reworking it to run really quickly on an entire array of timestamps.

Here is an example of the output of my range_diff function:

image

Using Irena's matlab scripts with ~similar settings, this is the result:

image

Focusing on specific time stamps, the two methods produce curves with similar shapes (though not the same due to slightly different window sizes - the matlab script requires some interpolation whereas ours doesnt) but completely different amplitudes. Here's ours:

image

and here's Irena's MATLAB results:

image

Additionally, our method creates this weird curved artifact for more distant times that does not show up in Irena's results. Seen here:

image

vs

image

Again, similar structure, especially near the basal reflector, but not quite there yet. This is also pre-phase unwrapping, which makes the results a lot odder.

Nevertheless, I think this is a start. Im now trying to just directly adapt the matlab method to python.

glugeorge commented 1 year ago

I've figured out the weird curving artifact. I was using wrong values for the speed of light in ice and the center wavelength in the phase2range function. I was tracking my code along with the MATLAB scripts pretty much step by step, and that was where the results diverged greatly. Using correct values fixes things.

So far it looks like the script in python (without wrapping) is pretty much spot on with Irena's scripts at its most basic form (which is the form that is the default setting). The one major difference is that the MATLAB scripts interpolate to ranges that are defined by a window size set by meters, whereas the python one is currently using indices for windowing. Now, the amplitudes of our data are a lot more similar, and the patterns are in close alignment:

image

This holds up pretty well even with the cumulative summation to get the total displacement after a period of time. After 100 time steps, it's still quite similar, with the exception of a few data points in the noisy sections:

image

Consequently, I'm pretty content with this being able to represent the raw xcor output from Irena's scripts - I think differences now are limited to what is essentially noise.

My next focus is making sure my phase wrapping is done similarly - in the MATLAB plots the phase wrapped results are in red, and they differ a lot from the current attempts to introduce phase wrapping in my scripts.

jkingslake commented 1 year ago

Is this the PR I should merge now to correct the data type bug? Or does it also contain the vertical velocity changes that are not ready yet?

glugeorge commented 1 year ago

I am happy to merge if you are, I think it is in an ok enough position with enough relevant changes to be incorporated in the main branch

glugeorge commented 1 year ago

Added in the README, I think merging is blocked on my end until you approve the changes?

jkingslake commented 1 year ago

thanks for all this! I have just merged it!