Re: Problem array subscripting

[Craig Markwardt <craigmnet(at)cow.physics.wisc.edu>]

GRI <dibeas@etsit.upm.es> writes:
> 
> I have a little problem I want to solve using array subscripting (if
> possible), but I haven't been able to do it without using a FOR
> statement:
> 
> I have a 3D array (dimX, dimY, n_images), which represents a time-series
> of images.
>
> I have a masking image (dimX, dimY), from wich I extract some points I
> am interested in (a ROI) with the command WHERE.
> 
> I could then change the one dimensional subscripts returned by WHERE
> into two vectors, that give me the X and Y coordinates of the points.
> 
> What I want to do is:
> 
> index= INDGEN(n_images)
> result= ARRAY[X, Y, index],

First of all, I would argue that using a FOR loop in this situation is
acceptable.  Why?  Because you are still vectorizing the majority of
the computation at the image level.  The extra computational overhead
of the FOR loop for a few images will be small.

That said, if you *really* want to vectorize, then I think you need to
take a different approach.  Array subscripting with index lists can
get tricky if you mix with other forms of subscripting at the same
time.  The WHERE command, operated on a two dimensional array, gives a
*one* dimensional list of array indices.  Therefore, it's usually best
to refer to your data the same way, ala REFORM.

If MASK is (NX x NY) and ARRAY is (NX x NY x N_IMAGES) then this
should work.

wh = where(MASK EQ 1, n_pix)  ;; One-dimensional array
array = reform(array, nx*ny,  n_images, /overwrite)
result = array(wh, *)
array = reform(array, nx, ny, n_images, /overwrite) ;; Revert to old dims.

It may look ugly, but it's pretty fast since the
REFORM(...,/OVERWRITE) operates on the data in place.  The dimensions
of RESULT are (N_PIX x N_IMAGES) where N_PIX is the number of pixels
selected by WHERE.

Good luck,

Craig


-- 
--------------------------------------------------------------------------
Craig B. Markwardt, Ph.D.         EMAIL:    craigmnet@cow.physics.wisc.edu
Astrophysics, IDL, Finance, Derivatives | Remove "net" for better response
--------------------------------------------------------------------------