| Re: phase unwrapping [message #26689 is a reply to message #20446] |
Tue, 18 September 2001 11:36   |
Clay Kirkendall
Messages: 5 Registered: August 2001
|
Junior Member |
|
|
Alexander Rauscher wrote:
> hi everybody
>
> i need a phase unwrapping algorithm (MRI, buti think it doenst really
> matter what type of data it is. does anyone have experience with
> this? it turnded out to be quite a tough problem...
> alex
Alex,
Here is a routine that I have been using for several years without
problems.
It defaults to 2pi phase resets but by using thresh and step it will
correct
for any step size. Depending on your data you may need to manually set
thresh and
step.
Good luck,
Clay
Function unwrap, x, thresh, step
;
; This routine removes resets from data vectors that have a
;modulus reset of size step. The data record is scanned for
;|resets| > |thresh| and adjusts the x by multiples of step
;
;
;
;MODIFICATIONS: 2/99 - added defaults for thresh and step
; 000831 - let x be an array
;
if n_params() eq 1 then $
Begin
thresh = !pi * 1.5
step = 2 * !pi
EndIf
sx=size(x)
if sx(0) EQ 0 then $
BEGIN
print, 'X must be a vector for this routine'
return, NaN
ENDIF
if sx[0] EQ 2 then $ ;Input is a two dimensional array
Begin
lx = length(x)
;get the derivative
dt=x(*, 0:lx-2) - x(*, 1:lx-1)
;search for steps greater then thresh
ddt=float(dt ge thresh) - (dt le (-thresh))
iddt=[ [replicate(0.0, sx[1])], [total(ddt, 2,
/CUMULATIVE)] ] * step
return, x + iddt
EndIf Else $
Begin
lx = max(sx(1:sx(0)), dim)
;get the derivative
dt=x(0:lx-2) - x(1:lx-1)
;search for steps greater then thresh
ddt=(dt ge thresh)*1. - (dt le (-thresh))
iddt=[0, total(ddt, /CUMULATIVE)] * step
return, x + iddt
EndElse
END
|
|
|
|