comp.lang.idl-pvwave archive: archive » Re: algorithm question. Can I get rid of the for loop?

Home » Public Forums » archive » Re: algorithm question. Can I get rid of the for loop?

Show: Today's Messages :: Show Polls :: Message Navigator
E-mail to friend

Re: algorithm question. Can I get rid of the for loop? [message #83715 is a reply to message #83706]

Wed, 27 March 2013 06:31

Søren Frimann
Messages: 13
Registered: November 2010

Junior Member

Den fredag den 22. marts 2013 14.25.44 UTC+1 skrev Heinz Stege:
>
> Before entering the loop get the indices of the lower and upper limits
>
> for all x values by use of the VALUE_LOCATE function. Then you can use
>
> this pre-calculated indices instead of the index array from the WHERE
>
> function.

That was a very nice hint indeed! Below an implementation where this has been done (as well as a few general updates). It's much faster than my older solution although, it retains the for loop

#########################################################

FUNCTION hampel, x, y, dx, THRESHOLD=threshold

Compile_Opt idl2

IF N_Elements(threshold) EQ 0 THEN threshold = 3

s0 = FltArr(N_Elements(y))
y0 = FltArr(N_Elements(y))
yy = y

lower_Boundary = Value_Locate(x,x-dx)+1 ; indices of lower boundaries
upper_Boundary = Value_Locate(x,x+dx) ; indices of upper boundaries

FOR i=0,N_Elements(y)-1 DO BEGIN
IF lower_Boundary[i] EQ upper_Boundary[i] THEN BEGIN
; only one point in gap
y0[i] = y[i]
s0[i] = !Value.F_NAN
ENDIF ELSE BEGIN
; Two or more points in gap
y_temp = y[lower_Boundary[i]:upper_Boundary[i]]
y0[i] = Median(y_temp) ; median filtering
s0[i] = 1.4826*Median(Abs(y_temp - y0[i])) ; estimating uncertainty
ENDELSE
ENDFOR

gp = Where(Abs(y - y0) LT threshold*s0) ;index of good points
ol = Where(Abs(y - y0) GE threshold*s0,n) ;index of outliers

yy[ol] = y0[ol] ; replace outliers

result = Create_Struct('y' ,yy, $
'sigma',s0, $
'gp' ,gp, $
'ol' ,ol, $
'n' ,n) ; number of outliers

RETURN, result

END

#########################################################

Den fredag den 22. marts 2013 18.29.31 UTC+1 skrev bobgst...@gmail.com:
> Using histogram and reverse indices would be much faster than looping and whereing. (i.e. get all of your "index" arrays in one call, rather than n_elements(y) calls of where).
>

I've looked into it, and I really don't see a way of using histogram, since it involves binning of the data, and my data aren't binned - rather they are subject to a moving window running smoothly over the data set.

Cheers,
Søren

Report message to a moderator

[Message index]

		Re: algorithm question. Can I get rid of the for loop? By: Jeremy Bailin on Wed, 27 March 2013 11:55
		Re: algorithm question. Can I get rid of the for loop? By: Jeremy Bailin on Wed, 27 March 2013 10:00
		Re: algorithm question. Can I get rid of the for loop? By: Jeremy Bailin on Wed, 27 March 2013 09:53
		Re: algorithm question. Can I get rid of the for loop? By: Søren Frimann on Wed, 27 March 2013 06:31
		Re: algorithm question. Can I get rid of the for loop? By: bobgstockwell on Fri, 22 March 2013 10:29
		Re: algorithm question. Can I get rid of the for loop? By: cgguido on Fri, 22 March 2013 10:01
		Re: algorithm question. Can I get rid of the for loop? By: Heinz Stege on Fri, 22 March 2013 06:25
		Re: algorithm question. Can I get rid of the for loop? By: Jeremy Bailin on Thu, 21 March 2013 20:39

Previous Topic:	Rotating Tickmark labels
Next Topic:	Re: Inscribed circle

-=] Back to Top [=-

[ Syndicate this forum (XML) ] [

] [

]

Current Time: Fri Oct 10 23:12:25 PDT 2025

Total time taken to generate the page: 0.00723 seconds