craigmnet@cow.physics.wisc.edu wrote:
> After looking at your problem it looks like the widths of the lines
> are as wide as your energy range. I guess this would be appropriate
> [snip..]
> so many exponentiations. You could have gotten a pretty big savings
> if the lines were narrow, and you could restrict the computation to a
> narrow region around the line center (say +/- 10 sigma with WHERE). I
> recommend that anyway to get rid of underflow errors.
> [snip..]
> If for example your sigma term were
> (((.12*sqrt(stick[0,i]/1000))/1.6651)*10) [ note the last factor is
> smaller ] then things start to look interesting.
That's because there's always a danger in taking code x and modifying it
into simpler example y! As you suggested, I actually DO have a sigma term
with a *10 instead of a *1000.
I did what you suggested with the +/- 10*sigma with where (code below),
and it drastically improved running time. For the 5000 element array, the
running time went from 62.3 seconds to 4.9 seconds! I plot the error
associated with the method, and I was actually able to go down to +/-
4*sigma before there was any noticable error.
Thanks!
Todd
--
pro test
convoluted = fltarr( 2, 2048 )
convoluted[0,*] = findgen( 2048 ) / 2047. * 4.
convoluted2 = convoluted
;; Let's fake a stick spectrum, we usually have at least this many elements
stick = abs(randomn(10, 2, 5000 ))
stick[0,*] = stick[0,*] * 4.
stick[1,*] = stick[1,*] * 1000.
time=systime(1)
for i=0L, n_elements( stick ) / 2 -1 do $
convoluted[1,*] = stick[1,i]*exp(-((convoluted[0,*] - stick[0,i])^2)/ $
(((.12*sqrt(stick[0,i]/1000))/1.6651)*10)^2) $
+ convoluted[1,*]
time2 = systime(1)
for i=0L, n_elements( stick ) / 2 -1 do begin
;; Find range
sigma = (((.12*sqrt(stick[0,i]/1000))/1.6651)*10)
range = where( convoluted2[0,*] ge (stick[0,i]-4*sigma) and $
convoluted2[0,*] le (stick[0,i]+4*sigma) )
if( range[0] ne -1 ) then $
convoluted2[1,range] = stick[1,i]*exp(-((convoluted2[0,range] - $
stick[0,i])^2)/sigma^2) $
+ convoluted2[1,range]
endfor
print, 'time old: ', time2-time
print, 'time new: ', systime(1)-time2
new = convoluted
new[1,*] = new[1,*] - convoluted2[1,*]
plot,new[0,*], new[1,*],yrange=[-10,10]
end ;; test
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
