four dimensional array

[Paola Madonna]

 

 

Hi,

we have problems with the function interpn in a c++ program when we have a four dimensional array; the results are different from those obtained with octave line command.

Using a three dimensional array defined  as  

 

NDArray lookupT(dim,0.0);     

for (octave_idx_type i=0;i<3;i++)

            for (octave_idx_type j=0;j<3;j++)

                for (octave_idx_type k=0;k<3;k++)

                    {

                        lookupT(i,j,k) = i+j+k;

                     }

the interpn works well as shown in the following code:

 

#include <octave/oct.h>

#include <octave/octave.h>

#include <octave/parse.h>

#include <iostream>

using namespace std;

 

int     main (int argc, char **argv)

{

 

    if (octave_main (argc, argv, 1))

    {

        dim_vector  dim;

        dim.resize(3);

        dim(0)=3;

        dim(1)=3;

        dim(2)=3;

 

        NDArray lookupT(dim,0.0);

        for (octave_idx_type i=0;i<3;i++)

            for (octave_idx_type j=0;j<3;j++)

                for (octave_idx_type k=0;k<3;k++)

                    {

                        lookupT(i,j,k) = i+j+k;

                     }

 

        dim_vector dimv;

        dimv.resize(2);

        dimv(0) = 3;

        dimv(1) = 1;

        NDArray nsv(dimv);

        NDArray cee(dimv);

        NDArray ro(dimv);

 

        for (octave_idx_type k=0;k<3;k++)

{

            nsv.elem(k)=k;

            cee.elem(k)=k;

            ro.elem(k)=k;

        }

 

           octave_value_list f_arg;

 

           f_arg(0) = octave_value(nsv);

           f_arg(1) = octave_value(cee);

           f_arg(2) = octave_value(ro);

           f_arg(3) = octave_value(lookupT);

           f_arg(4) = 1.1;

           f_arg(5) = 1.5;

           f_arg(6) = 2.0;

           octave_value_list retval = feval("interpn",f_arg);

           if (retval.length()>0)

           {

               Matrix res(retval(0).matrix_value());

               cout<<"RES "<<res<<endl;

           }

           else

               error ("feval failed");

    }

    else

        error ("Octave interpreter initialization failed");

    return error_state ? 1 : 0;

}

 

Otherwise using a three dimensional array defined as

      

        Array<idx_vector> idxt(3);

        idxt(0)=idx_vector(Range(0.,2));

        idxt(1)=idx_vector(Range(0.,2));

        idxt(2)=idx_vector(Range(0.,2));

for (octave_idx_type i=0;i<3;i++)

            for (octave_idx_type j=0;j<3;j++)

                for (octave_idx_type k=0;k<3;k++)

                    {

                        lookupT(idxt) = i+j+k;

                    }

 

The interpn doesn’t work well as shown in this code:

 

#include <octave/oct.h>

#include <octave/octave.h>

#include <octave/parse.h>

#include <iostream>

#include "time.h"

using namespace std;

 

int     main (int argc, char **argv)

{

 

    if (octave_main (argc, argv, 1))

    {

        dim_vector  dim;

        dim.resize(3);

        dim(0)=3;

        dim(1)=3;

        dim(2)=3;

 

        NDArray lookupT(dim,0.0);

        Array<idx_vector> idxt(3);

 

        idxt(0)=idx_vector(Range(0.,2));

        idxt(1)=idx_vector(Range(0.,2));

        idxt(2)=idx_vector(Range(0.,2));

 

        for (octave_idx_type i=0;i<3;i++)

            for (octave_idx_type j=0;j<3;j++)

                for (octave_idx_type k=0;k<3;k++)                   

                    {

                        lookupT(idxt) = i+j+k;

                    }

 

        dim_vector dimv;

        dimv.resize(2);

        dimv(0) = 3;

        dimv(1) = 1;

        NDArray nsv(dimv);

        NDArray cee(dimv);

        NDArray ro(dimv);

for (octave_idx_type k=0;k<3;k++)

        {

            nsv.elem(k)=k;

            cee.elem(k)=k;

            ro.elem(k)=k;

        }

 

           octave_value_list f_arg;

 

           f_arg(0) = octave_value(nsv);

           f_arg(1) = octave_value(cee);

           f_arg(2) = octave_value(ro);

           f_arg(3) = octave_value(lookupT);

           f_arg(4) = 1.1;

           f_arg(5) = 1.5;

           f_arg(6) = 2.0;

           octave_value_list retval = feval("interpn",f_arg);

           if (retval.length()>0)

           {

               Matrix res(retval(0).matrix_value());

               cout<<"RES "<<res<<endl;

           }

           else

               error ("feval failed");

    }

    else

        error ("Octave interpreter initialization failed");

    return error_state ? 1 : 0;

}

 

Why does not interpn work well in the second code? Is it a possible known bug?

We need a four dimensional array and in this case we must define the four dimensional array as defined in the second code.

Could you give us an explanation?

 

Thanks in advance

Paola Madonna