How to use modified Shepard method

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Hi

It is my first using of AlgLib. So sorry if I'll wrong.

I have nonuniform mesh with data in some point and I need to recieve approximated value in other point of mesh. And if I understand right I need to use Shepard method ("idwbuildmodifiedshepard") for my purpose.

I have 2D arrary m_Array with double values

Code:

double m_Array[iArrayWidth][iArrayHeight];
//.....
// filling m_Array with value
//....

Then I've created:

Code:

ap::real_2d_array xy;
xy.setlength(iArrayWidth, iArrayHeight);

for ( int i = 0; i < iArrayHeight; i++ )
{
   for ( int j = 0; j < iArrayWidth; j++ )
   {
      xy(i,j) = m_Array[i][j];
   }
}

// pass to the idwbuildmodifiedshepard

Is it true when I fill xy array with my values?

Thanks a lot.
Best regards, Aleksey.

Sergey.Bochkanov · **Joined:** Fri May 07, 2010 7:06 am **Posts:** 906

You should store your points in rows, one row = one point. First iArrayWidth-1 columns store coordinates, last column stores function value.

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Thanks for reply! ;)
Did you mean that should I write something like that

Code:

ap::real_2d_array xy;
xy.setlength(m_iArrayHeight, m_iArrayWidth);

for ( int i = 0; i < m_iArrayHeight; i++ )
{
   for ( int j = 0; j < m_iArrayWidth; j++ )
   {
      xy(i*j,0) = XXX;
      xy(i*j,1) = m_ppArray[i][j].dValue;
   }
}

And if it true, what is XXX ?

Thanks
Best regards, Aleksey.

Sergey.Bochkanov · **Joined:** Fri May 07, 2010 7:06 am **Posts:** 906

It should be something like

Code:

ap::real_2d_array xy;
xy.setlength(NumberOfPoints, NumberOfCoordinates+1); // NumberOfCoordinates=3, for example

for ( int i = 0; i < NumberOfPoints; i++ )
{
    xy(i,0) = x[i]; // x-value
    xy(i,1) = y[i]; // y-value
    xy(i,2) = z[i]; // z-value
    xy(i,3) = f[i]; // function value
}
}

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Hi, Sergey. Thank you very much! :) It is now clearly for me.

Now I fill my array with following way

Code:

xy.setlength(m_iArrayHeight * m_iArrayWidth, 3);

int k = 0;
for ( int i = 0; i < m_iArrayHeight; i++ )
{
   for ( int j = 0; j < m_iArrayWidth; j++ )
   {
      xy(k, 0) = m_ppArray[i][j].point.x;
      xy(k, 1) = m_ppArray[i][j].point.y;
      xy(k, 2) = m_ppArray[i][j].dValue;
      k++;
   }
}

Now I have another question. Then I pass "xy" to "idwbuildmodifiedshepard" method, which returns IDW interpolant - Z. And if I see right, Z will contain approximated values. So which field of Z contains these values? I've tried following code:

Code:

k = 0;
for ( int i = 0; i < m_iArrayHeight; i++ )
{
   for ( int j = 0; j < m_iArrayWidth; j++ )
   {
      double val = Z.q(k++, 2);
   }
}

But values of "val" are the same as in m_ppArray[i][j].dValue, i.e. it no changed. Maybe I've made something wrong again?

Thanks a lot!
Best regards, Aleksey.

Sergey.Bochkanov · **Joined:** Fri May 07, 2010 7:06 am **Posts:** 906

You should use idwcalc() function to calculate the value of the interpolant at the points you need. Never try to work with its fields directly.

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Hi Sergey. Thank you again! :)

To calculate the value of the interpolant I use following code:

Code:

ap::real_1d_array x;
idwinterpolant z1;
int n = m_iArrayWidth * m_iArrayHeight;

idwbuildmodifiedshepard(xy, n, 2, 1, 15, 25, z1);

x.setlength(nx);

for( int i = 0; i < n; i++ )
{
   x(0) = xy(i, 2);
   //ap::vmove(&x(0), 1, &xy(i, 2), 1, ap::vlen(0,nx-1));
   double val = idwcalc(z1, x);
   idwerrors = ap::fp_neq(val, xy(i,nx));
}

In my case method idwcalc almost always returns 0 and sometimes -1.7054547436862421e+073. So I think I've made something wrong again. Could you help me once again please?

Thanks
Best regards, Aleksey.

Sergey.Bochkanov · **Joined:** Fri May 07, 2010 7:06 am **Posts:** 906

Can you post here (as zip-file) your program with all its data? It is hard to tell what's wrong without looking into it. I think that IDW is too robust to produce such large deviations (1E+73), but I have to look at it myself.

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Hi, Sergey

I've added to archive "signal.txt" file with my data and cpp file with one method which reads signal.txt file and calcs interpolant.

Thank for your attention paid!
Best regards, Aleksey.

AlekseyL · **Joined:** Mon Sep 13, 2010 12:58 am **Posts:** 14

Hey Sergey.

Did you see my example code? Is it wrong? Or may be did you want something else?

Thanks!
Best regards, Aleksey.

forum.alglib.net

Forum rules

How to use modified Shepard method

Who is online