Hi,
I'm vectorizing a procedure in order to delete three for loops and
increase the speed.
I have to calculating the interpolation on a regular grid.
newy=interpol(y,x,newx,/spline)
 in the loop manner  for each index of the loop  there are  three
different arrays .
in the vectorize manner I have three 3D matrices, But I do not found a
way to calculate the interpolation matrix.
Any one have suggestions?

Thanks

Romolo

In article 
<041dea05-b2af-435a-b9e9-7ca360287c69@p35g2000yqh.googlegroups.com>,
 Romolo Politi <romolo.politi@gmail.com> wrote:

> Hi,
> I'm vectorizing a procedure in order to delete three for loops and
> increase the speed.
> I have to calculating the interpolation on a regular grid.
> newy=interpol(y,x,newx,/spline)
>  in the loop manner  for each index of the loop  there are  three
> different arrays .
> in the vectorize manner I have three 3D matrices, But I do not found a
> way to calculate the interpolation matrix.
> Any one have suggestions?
> 
> Thanks
> 
> Romolo

I'm afraid that you have not clearly explained your problem.  Are you
trying to do 3-D cubic spline interpolation?  As far as I know,
standard IDL does not contain a multi-dimensional spline interpolation
procedure, only the 1-D functionality available in INTERPOL, SPLINE, and
SPL_INTERP.  You can do vectorized multi-dimensional *linear* interpolation
by using INTERPOLATE.

The Advanced Math and Statistics package (extra cost) includes the IMSL
spline interpolation functions.  You can view the help files and
see if that is what you need.

   http://127.0.0.1:60523/help/index.jsp

Ken Bowman

On 10 Dic, 16:25, "Kenneth P. Bowman" <k-bow...@null.edu> wrote:
> In article
> <041dea05-b2af-435a-b9e9-7ca360287...@p35g2000yqh.googlegroups.com>,
> =A0Romolo Politi <romolo.pol...@gmail.com> wrote:
>
> > Hi,
> > I'm vectorizing a procedure in order to delete three for loops and
> > increase the speed.
> > I have to calculating the interpolation on a regular grid.
> > newy=3Dinterpol(y,x,newx,/spline)
> > =A0in the loop manner =A0for each index of the loop =A0there are =A0thr=
ee
> > different arrays .
> > in the vectorize manner I have three 3D matrices, But I do not found a
> > way to calculate the interpolation matrix.
> > Any one have suggestions?
>
> > Thanks
>
> > Romolo
>
> I'm afraid that you have not clearly explained your problem. =A0Are you
> trying to do 3-D cubic spline interpolation? =A0As far as I know,
> standard IDL does not contain a multi-dimensional spline interpolation
> procedure, only the 1-D functionality available in INTERPOL, SPLINE, and
> SPL_INTERP. =A0You can do vectorized multi-dimensional *linear* interpola=
tion
> by using INTERPOLATE.
>
> The Advanced Math and Statistics package (extra cost) includes the IMSL
> spline interpolation functions. =A0You can view the help files and
> see if that is what you need.
>
> =A0 =A0http://127.0.0.1:60523/help/index.jsp
>
> Ken Bowman

Sorry for the not clarity.
my situation is
x=3Dfltarr(3200,201,201)
y=3Dfltarr(3200,201,201)
newx=3Dfltarr(799,201,201)
newy=3Dfltarr(799,201,201)
for i=3D0,200 do begin
  for j=3D0,200 do begin
    newy[*,i,j]=3Dinterpol(y[*,i,j],x[*,i,j],newx[*,i,j],/SPLine)
  endfor
endfor

I like vectorize this loop.

Thanks, and sorry again

Romolo Politi

In article 
<b599bcdf-d42e-45d8-b643-7ea4620fd7ed@k19g2000yqc.googlegroups.com>,
 Romolo Politi <romolo.politi@gmail.com> wrote:

> Sorry for the not clarity.
> my situation is
> x=fltarr(3200,201,201)
> y=fltarr(3200,201,201)
> newx=fltarr(799,201,201)
> newy=fltarr(799,201,201)
> for i=0,200 do begin
>   for j=0,200 do begin
>     newy[*,i,j]=interpol(y[*,i,j],x[*,i,j],newx[*,i,j],/SPLine)
>   endfor
> endfor

You are doing ~40,000 separate spline fits.  I don't think there is
any way to vectorize this, as each fit is an independent problem.
Spline fitting is a coupled problem, so you cannot, for example, simply
string your segments together into one long fit.

If the newx[*,i,j] are the same for all i and j, you can probably get
a significant speed-up by calling SPL_INIT once and then calling
SPL_INTERP for each i and j.  Since you are storing newx separately
for each i and j, I gather this is not the case.

Your program is accessing memory efficiently, but you will get a small
speed up by saying

   newy[0,i,j]=interpol(y[*,i,j],x[*,i,j],newx[*,i,j],/SPLine)

You might also try the other built in spline functions:  SPLINE,
SPLINE_P, SPL_INIT and SPL_INTERP to see whether one is faster.

Finally, I notice that you are downsampling your data substantially.
Are you sure you really need to use a spline fit to do that?

Ken Bowman

On Dec 10, 3:00=A0pm, "Kenneth P. Bowman" <k-bow...@null.edu> wrote:
> In article
> <b599bcdf-d42e-45d8-b643-7ea4620fd...@k19g2000yqc.googlegroups.com>,
> =A0Romolo Politi <romolo.pol...@gmail.com> wrote:
>
> > Sorry for the not clarity.
> > my situation is
> > x=3Dfltarr(3200,201,201)
> > y=3Dfltarr(3200,201,201)
> > newx=3Dfltarr(799,201,201)
> > newy=3Dfltarr(799,201,201)
> > for i=3D0,200 do begin
> > =A0 for j=3D0,200 do begin
> > =A0 =A0 newy[*,i,j]=3Dinterpol(y[*,i,j],x[*,i,j],newx[*,i,j],/SPLine)
> > =A0 endfor
> > endfor
>
> You are doing ~40,000 separate spline fits. =A0I don't think there is
> any way to vectorize this, as each fit is an independent problem.
> Spline fitting is a coupled problem, so you cannot, for example, simply
> string your segments together into one long fit.
>
> If the newx[*,i,j] are the same for all i and j, you can probably get
> a significant speed-up by calling SPL_INIT once and then calling
> SPL_INTERP for each i and j. =A0Since you are storing newx separately
> for each i and j, I gather this is not the case.
>
> Your program is accessing memory efficiently, but you will get a small
> speed up by saying
>
> =A0 =A0newy[0,i,j]=3Dinterpol(y[*,i,j],x[*,i,j],newx[*,i,j],/SPLine)
>
> You might also try the other built in spline functions: =A0SPLINE,
> SPLINE_P, SPL_INIT and SPL_INTERP to see whether one is faster.
>
> Finally, I notice that you are downsampling your data substantially.
> Are you sure you really need to use a spline fit to do that?
>
> Ken Bowman

You might be able to fake something in the "stringing the segments
together" vein by putting some padding in between the segments that
does a linear interpolation between the endpoint of one segment and
the beginning of the next, enough padding to make the last old segment
and first new segment independent. Of course, you'll get a different
answer for those end segments because of that padding than if they
didn't exist, so that may not work depending on how much those points
are offset, but it should at least be better than sticking them next
to each other. What boundary conditions does INTERPOL use with spline
fits? Maybe you can design padding that effectively gives you the same
boundary conditions?

-Jeremy.