comp.lang.idl-pvwave archive
Messages from Usenet group comp.lang.idl-pvwave, compiled by Paulo Penteado

Home » Public Forums » archive » dimensional headache
Show: Today's Messages :: Show Polls :: Message Navigator
E-mail to friend 
Switch to threaded view of this topic Create a new topic Submit Reply
dimensional headache [message #43322] Fri, 01 April 2005 08:49
Margrethe is currently offline  Margrethe
Messages: 4
Registered: March 2005
Junior Member
I'm developing a two-dimensional model and to save computational time
I've tried to avoid for-loops whereever I could and used matrices
instead. I thought I had it all figured out, but I'm doing something
wrong...I'd be grateful if someone could tell me what's wrong with the
following:

The model lies on a grid xgrid=findgen(160), ygrid=findgen(160). I want
to calculate the numerical value of the model at the coordinate
xmod,ymod

Here's my old code using for loops

for i=0,159 do begin
for j=0,159 do begin
dx_sqr = ( xmod - xgrid (i) )^2.
dy_sqr = ( ymod - ygrid (j) )^2.
d = sqrt ( dx_sqr + dy_sqr )
; don't take the wings of the psf into account
indx = where (d lt (3. * (fwhm/2.35) ) ,ct)

if (ct gt 0) then begin
vel_p = ( vel_sys + v_rad(i,j) )
psf = gaussfunc(st_dev,xgrid(i),ygrid(j),xmod,ymod )
lineprof(i,j) = lineprof(i,j) + $
total ( psf * image(i,j) * $
exp( (vel - vel_p(i,j))^2./(-2.*v2)) )
endif

endfor
endfor

The above example works fine, but of course takes a long time to
compute.
Here's my attempt with matrices, which doesn't give the correct output.


dx_sqr = (xmod - findgen(160))^2.
dy_sqr = (ymod - findgen(160))^2.

; expand one dimension
dx_sqr_expand = rebin(dx_sqr,160,160)
dy_sqr_expand = transpose(rebin(dy_sqr,160,160))

d = sqrt(dx_sqr_expand + dy_sqr_expand)

; don't take the wings of the psf into account

WhereToMulti,d,where(d lt (3.*(fwhm/2.35)),ct),colindx,rowindx

vel_p = (vel_sys + v_rad(colindx,rowindx))

psf = gaussfunc(st_dev,xgrid(colindx),ygrid(rowindx),xmod,ymod)

lineprof(k,l) = lineprof(k,l) + $
total( psf * image(colindx,rowindx) * $
exp( (vel - vel_p(colindx,rowindx))^2./(-2.*v2)) )
  Switch to threaded view of this topic Create a new topic Submit Reply
Previous Topic: Re: IDLWAVE Manuals: Your Support Needed
Next Topic: Last rites declaration of Ioannes Paulus PP. II (Karol Wojtyla)

-=] Back to Top [=-
[ Syndicate this forum (XML) ] [ RSS ] [ PDF ]

Current Time: Wed Oct 08 15:34:30 PDT 2025

Total time taken to generate the page: 0.00441 seconds