| Re: Programming style [message #11740 is a reply to message #3375] |
Fri, 15 May 1998 00:00  |
J.D. Smith
Messages: 214
Registered: August 1996
|
Senior Member |
|
|
Patrick V. Ford wrote:
>
> I have a general style/algorithm question.
>
> I want to plot in a 3-D array an ellipsoid within an in ellipsoid where the
> voxels between the boundaries are non-zero and else where zero.
>
> The general function is (x/a)^2 + (y/b)^2 + (z/c)^2 = 1.0
>
> I have already done this using for loops and conditional statements but it
> occurred to me that there may be some IDL matrix-boolean logic combination
> that could accomplish this in a faster and more 'elegant' fashion.
>
> I am now open for suggestions?
>
> Thanks in advance.
>
> Patrick Ford, MD
> Department of Radiology
> Baylor College of Medicine
> pford@bcm.tmc.edu
How about:
recast as:
Eq. 1: z^2+(x/e11)^2+(y/e21)^2 z^2=c1^2 (outer ellipse)
Eq. 2: z^2+(x/e12)^2+(y/e22)^2 z^2=c2^2 (inner ellipse)
let your array be nx by ny by nz.
Then:
z=findgen(nx*ny*nz)
x2=((z mod (nx*ny)) mod nx-nx/2)^2
y2=((z mod (nx*ny))/nx-ny/2)^2
z2=((temporary(z)/(nx*ny)) - nz/2)^2
array=bytarr(nx,ny,nz)
array[where(z2+x2/e11^2+y2/e21^2 le c1^2 AND z2+x2/e12^2+y2/e22^2 ge
c2^2)]=1b
Note the two ellipses are centered on the midpoint of the array and are
concentric. This can be modified by changing the subtracted value in
each of x2,y2,z2. Definitely faster than loops. Elegance is in the eye
of the beholder, though.
NB: The x,y, and z index vectors must be floats, since for 3-d
data,indices get large pretty quick. E.g. 100x100x100 would choke with
longs (since 100^3^2=10^12=2^39.86!). This introduces some "fuzziness"
at the boundaries due to roundoff. You can throw in a floor() statement
to eliminate this if you really want.
JD
--
J.D. Smith |*| WORK: (607) 255-5842
Cornell University Dept. of Astronomy |*| (607) 255-4083
206 Space Sciences Bldg. |*| FAX: (607) 255-5875
Ithaca, NY 14853 |*|
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
