cgrima / subradar

Tools for subsurface radar processing
MIT License
3 stars 3 forks source link

gcc iteration #6

Open gregoryng opened 1 year ago

gregoryng commented 1 year ago

The code around the gcc iteration is a little bit peculiar.

https://github.com/cgrima/subradar/blob/ca50a1cb5f4a2196d202d921a7b69eb5b873011a/subradar/surface.py#L172-L177

The main loop appears to iterate i through the interval [0, ysize-2], so that you calculate the correlation of each trace rdg[i] with the following trace rdg[i+1], obviously stopping before you go past the end of the radargram.

Then, for the final case, you use the loop variable i, and it seems like you intend to use the last trace and the one before it, and do the same thing.

https://github.com/cgrima/subradar/blob/ca50a1cb5f4a2196d202d921a7b69eb5b873011a/subradar/surface.py#L183-L184

Would it be safe to say that this is another way of expressing your intent?

 _ = utils.gcc(rdg[:, -2], rdg[:, -1], **kwargs)

The problem is that unlike in C/C++, the loop variable i doesn't advance past the end of the range to ysize. It remains at the final value in the loop, i.e., yn[-2]. So that makes this code not quite do what you want.

For example, if you have a radargram with 100 traces, then the above loop will iterate i through 99 values -- [0...98]. The final value will be 98, not 99.

However there is another optimization. If your intent is to use the last trace and the next-to-last trace, then you actually do not need to recalculate utils.gcc again, because you already did that calculation in the last iteration of the for loop, and placed it into the next-to-last elements of your output vectors. You should be able to go like this:

    # Last record
    # _ = utils.gcc(rdg[:, i], rdg[:, i-1], **kwargs)
    tau[-1] = tau[-2] 
    val[-1] = val[-2]
    cc[:,-1] = cc[:, -2]
    #ch[:,-1] = ch[:, -2]

If you agree with this, I can roll this into some code I am working on and submit a pull request. There was just enough uncertain here that I wanted to make sure I understood your idea and explained what I was thinking.

Oh one other question -- did I get the axes right? Is axis 0 fast time and axis 1 slow time? Or the other way around? I can add notes in the comments to clarify that.

cgrima commented 1 year ago

Thanks for reviewing that. The cross-correlation (gcc) technique is still a WIP/exploratory approach that I never use in production for detecting the surface for now. I will fix those issues once I have some time to focus back on this technique.

gregoryng commented 1 year ago

The code has some similarities to what Kirk was doing for interferometry in SDS -- though Kirk also wanted to allow for sub-sample cross correlation and roll.