ダイアログベースの3次Spline補間のプログラム

このQ&Aのポイント
  • Visual C++ 2008 Express Editionを使用して、ダイアログベースの3次Spline補間のプログラムを作成しています。
  • 各節点間を等分割することはできましたが、各区間を同一刻み幅で算出する方法がわかりません。
  • プログラムの一部を示していますが、dnの配列に入ってくる分割数の値に違いがあることがわかります。
回答を見る
  • ベストアンサー

ダイアログベースの3次Spline補間のプログラム

現在、Visual C++ 2008 Express Editionを使用して下に示すようなダイアログベースの3次Spline補間のプログラムを作成しています。 各節点間を等分割するほうはできたのですが、各区間を同一刻み幅で算出するほうが、なかなかできません。 どこをどうしたらいいかがさっぱりです。 以下にプログラムの一部(各区間を同一刻み幅で計算する部分)を示しますので、教えてください。 わからない点は、dnの配列に入ってくる分割数が合っているときと少ないときがあることです。プログラムの後に使用した入力データと、dnの中身を参考までに掲載します。 ============▼プログラムの一部=========== /* n:全点数(int), num:全節点数(int),xn,yn:節点(double型のarray),dX:全区間(xn[0]~xn[num-1])の刻み幅(double) * dn[]:各区間の分割数(int),xx,yy:Spline補間による点(節点を含む,double),x計算された補間点のx座標(double) Cub(double):3乗する関数,Sqr(double):2乗する関数,a[],b[],c[],d[]:3次Spline補間に用いる3次関数の各係数*/ //nの初期化 n =0; //テキストボックスから刻み幅を取得 dX = double::Parse(DX->Text::get()); //各区間の分割数を計算し、全点数を求める for(i=0;i<num-1;i++){ dn[i]=(int)((xn[i+1]-xn[i])/dX); n +=dn[i]; } //xx,yyを全点数+1(最後の節点を含むために+1をしている)で領域確保 xx=gcnew array<double>(n+1); yy=gcnew array<double>(n+1); //ループインデックスを初期化 i=0,j=0; //補間点のx座標を求めながらy座標を算出 for(x=xn[0];x<=xn[num-1];x+=dX){ if(!((xn[j]<=x)&&(x<xn[j+1]))){ j++; } if(i>n+1){ break; } xx[i]=x; yy[i]=a[j]*Cub(xx[i]-xn[j])+b[j]*Sqr(xx[i]-xn[j])+c[j]*(xx[i]-xn[j])+d[j]; i++; } =======▼入力データ(csvファイル)=========== xn yn 0.904 7.54 1.002 13.86 1.104 23.42 1.207 36.61 1.305 52.91 1.403 73.4 1.505 97 1.603 123 1.706 151.3 1.808 182.6 1.906 215.3 2.009 249.2 2.107 284.6 2.209 321.6 2.312 359 2.41 398.1 2.571 438.7 2.669 504.9 2.772 548.3 2.869 593 2.972 639 3.075 687 3.172 736 3.275 787 3.378 840 3.476 895 3.578 951 3.676 1011 3.779 1073 3.876 1136 3.979 1204 4.077 1276 4.18 1353 4.277 1440 4.38 1543 4.483 1762 4.585 2082 4.683 3445 4.786 5598 4.888 7120 4.991 8080 =======▼dnの中身=========== Excelでの計算値 プログラムでの計算値 --------------------------------------- 98 97 102 102 103 102 98 97 98 98 102 101 98 98 103 102 102 102 98 97 103 102 98 98 102 101 103 102 98 98 161 161 98 97 103 102 97 97 103 102 103 103 97 96 103 102 103 103 98 97 102 101 98 98 103 102 97 96 103 103 98 97 103 102 97 97 103 102 103 102 102 102 98 97 103 102 102 102 103 102 上の表の右側はダイアログ上で、「区間間隔」をクリック(各ラジオボタンをクリックすると右横のテキストボックスが入力できるようになります。)してテキストボックスに0.001を入れた時の結果です。

質問者が選んだベストアンサー

  • ベストアンサー
  • redfox63
  • ベストアンサー率71% (1325/1856)
回答No.1

intに変換する際の浮動小数点の誤差切捨てが発生しているのでしょう for(i=0;i<num-1;i++){   dn[i]=(int)Math::Round( (xn[i+1]-xn[i])/dX), 0 );   n +=dn[i]; } といった具合で丸めてみましょう

yf491224
質問者

お礼

早速の回答ありがとうございます。 試してみます^^

yf491224
質問者

補足

あれから、ソ-スコードをいわれたとおり、修正したら無事できました。  小数点のところで、誤差が出ていたなんて思いもよらなかったです。 本当にありがとうございました。

関連するQ&A

  • ラグランジュの補間法のCプログラム

    昨日学校でラグランジュの補間法の問題をC言語のプログラムで解けという課題が出されました しかし、友達と相談してもよくわかりませんでした 課題は以下の問題です sin関数6点、(0.92+0.01x)、x=0,1,2,3,4,5を求めて、ラグランジュの方法でsin(0.923)を計算せよ ちなみに答えは、0.79742です 先生からサンプルのプログラムをもらいました 以下のサンプルプログラムを参考にして解いてくださいと言われたのですが、どうしても解けません すいませんが分かる方、よろしくお願いします #include <stdio.h> #include <math.h> #define N 6 //データ数 double x[N]={ 0.0,1.0,2.0,3.0,3.1,5.0}; //X座標 double y[N]={0.0,1.1,2.5,4.0,4.1,5.0}; //Y座標 double lagrange( double); int main() { double xx,yy; //補間計算 printf("XX\t\tYY\n"); for( xx=0.0; xx<=5.0; xx+=.2 ) { yy = lagrange( xx); printf("%8.2lf\t%8.2lf\n", xx, yy ); } return 0; } //補間サブルーチン double lagrange( double xx ) { double z[N]; double yy=0.0; int i,j; for( i=0; i<N; i++ ) { z[i] = 1.0; //係数計算 for( j=0; j<N; j++ ) if( i!=j ) z[i]*=(xx-x[j])/(x[i]-x[j]); //補間値計算 yy+=z[i]*y[i]; } return yy; } 上記はあくまでサンプルプログラムなので、中に入っている数値は適当です よろしくお願いします

  • スプライン補間

    x=[-1,0,1,2],y=[0,1,0,0]のデータで 区間x=0~1 をスプライン補間で計算させています。 MuPAD でcubicSplineを用いた場合と C言語によるアルゴリズム辞典から作ったソフトでは計算結果が 微妙に異なります。 どちらが3次スプライン補間として正しいのかお教え願えないでしょうか? あるいはどちらも正しいとして、スプラインの種別が違うのでしょうか? 非常に漠然としていますが、よろしくお願いします。 「自分のツールだとこういう結果だった」というようなアドバイスでも大歓迎です。 X      MuPAD        C言語によるアルゴ 0      1            1 0.125  0.922851563  0.9488281 0.25   0.8203125    0.853125 0.375  0.698242188  0.7246094 0.5    0.5625       0.575 0.625  0.418945313  0.4160156 0.75   0.2734375    0.259375 0.875  0.131835938  0.1167969

  • 線形補間

    線形補間での求め方 問題文: 1.数値を読み込む 2.xを読み込む 3.x<x1 または x>xnならエラー 3.x1<x<i+1 となるiを見つける 4.補間公式でyを求める 5.結果をプリントする #include<stdio.h> float hokan(void); int xn[] = {0,5,10,15,20,25,30,35,40,45,50,55,60,65,70,75}; double yn[] = {0.000,0.087,0.173,0.258,0.342,0.422,0.500,0.573, 0.624,0.707,0.766,0.819,0.866,0.906,0.939,0.965}; int x=5; void main() { printf("y = %lf\n", hokan()); } float hokan(void) { int i; double y; if(x <0|| 75<= x){ printf("エラーです\n"); } else{ for(i=0; xn[i] < x; i++) y = (yn[i+1] - yn[i]) * (x - xn[i]) / (xn[i+1] - xn[i]) + yn[i]; return y; } } 数表を最初にxnとynで表記しています。 このプログラムで線形補間が行われてないそうなのですが・・ 何か誤りがある様でしたらどなたか教えてください。

  • エクセルにて、3次スプラインのマクロを組みました。

    エクセルにて、3次スプラインのマクロを組みました。 結果がところどころしか表示されません。 何度見直しても自分では間違いがわかりません。 どなたか、間違い箇所がわかりませんでしょうか。 A列には5行目から、0-720までの721個の少数値 0 1.005586592 2.011173184 B列には5行目から、721個の少数値 5.83 5.69 5.66 D列には5行目から、0-720まで、1ずつ増加の整数値 0 1 2 E列には補間値をマクロから代入 以下マクロです。 お手数をおかけしますが、目を通して指摘していただけませんでしょうか。 Sub 補間_3次スプライン() Dim i As Long Dim h(1000) As Double Dim dif1(1000) As Double Dim dif2(1000) As Double Dim data_count As Long Dim lngDataEnd As Long Dim dataX(1000) As Double Dim dataY(1000) As Double Dim x, y, yy0, yy1, yy2, yy3 As Double lngDataEnd = ThisWorkbook.Sheets("3次スプライン").Range("A65536").End(xlUp).Row data_count = lngDataEnd - 5 For i = 1 To data_count dataX(i) = ThisWorkbook.Sheets("3次スプライン").Cells(i + 4, 1) dataY(i) = ThisWorkbook.Sheets("3次スプライン").Cells(i + 4, 2) Next i h(0) = 0 dif2(0) = 0 On Error Resume Next For i = 1 To data_count h(i) = dataX(i) - dataX(i - 1) '//間隔を計算 dif1(i) = (dataY(i) - dataY(i - 1)) / h(i) '//一次微分を計算 Next i For i = 1 To data_count '二次微分を計算 dif2(i) = (dif1(i + 1) - dif1(i)) / (dataX(i + 1) - dataX(i - 1)) Next i i = 1 For x = 0 To 720 If x < dataX(i) Then yy0 = dif2(i - 1) / (6 * h(i)) * (dataX(i) - x) * (dataX(i) - x) * (dataX(i) - x) '第1項 yy1 = dif2(i) / (6 * h(i)) * (x - dataX(i - 1)) * (x - dataX(i - 1)) * (x - dataX(i - 1)) '第2項 yy2 = (dataY(i - 1) / h(i) - h(i) * dif2(i - 1) / 6) * (dataX(i) - x) '第3項 yy3 = (dataY(i) / h(i) - h(i) * dif2(i) / 6) * (x - dataX(i - 1)) '第4項 y = yy0 + yy1 + yy2 + yy3 ThisWorkbook.Sheets("3次スプライン").Cells(x + 4, 5) = y Else: i = i + 1 End If Next x End Sub

  • C言語のプログラムで質問です。

    C言語のプログラムで質問です。 下のプログラム(最小二乗法の計算)を実行したところ -1.#IND00 というエラーが出てしまいます。 どこを直せばいいのでしょうか、教えてください。 #include <stdio.h> #include <math.h> /* gauss33.c */ #define N 3 main(){ double A[N][N],Aa[N][N]; double b[N],x[N], bb[N], e[N]; int n=N; int i, j, k; double akk, aik, s; double y[N]; double xx,yy; for(i=0;i<n;i++){ /*変数の初期化*/ x[i]=y[i]=0; for(j=0;j<n;j++) A[i][j]=0; } for(i=0;i<5;i++){ /*データ点は5点*/ printf("\n(x,y)="); scanf("%lf,%lf",&xx,&yy); A[0][0]+=xx*xx*xx*xx; /*Σx^4*/ A[0][1]+=xx*xx*xx; /*Σx^3*/ A[0][2]+=xx*xx; /*Σx^2*/ A[0][1]=A[1][0]; A[0][2]=A[1][1]=A[2][0]; A[1][2]+=xx; /*Σx*/ A[1][2]=A[2][1]; A[2][2]=n; y[0]+=xx*xx*yy; /*Σx^2y*/ y[1]+=xx*yy; /*Σxy*/ y[2]+=yy; /*Σy*/ } /* save original coefficients */ for(i=0; i<n; i++){ for(j=0; j<n; j++){ Aa[i][j]=A[i][j]; } bb[i]=b[i]; } /* forward operation */ for(k=0; k<n-1; k++){ akk=1/A[k][k]; for (i=k+1; i<n; i++){ aik=-A[i][k]*akk; for (j=k+1; j<n; j++){ A[i][j]+=aik*A[k][j]; } b[i]+=aik*b[k]; } for(j=k+1; j<n; j++){ A[k][j]*=akk; } b[k]*=akk; } /* backward operation */ x[n-1]=b[n-1]/A[n-1][n-1]; for(k=n-2; k>=0; k--){ s=0.0; for (j=k+1; j<n; j++){ s+=A[k][j]*x[j]; } x[k]=b[k]-s; } /* chek */ for(i=0; i<n; i++){ s=0.0; for(j=0; j<n; j++){ s+=Aa[i][j]*x[j];} e[i]=s-bb[i]; printf("\nx(%d)=%f error=%f\n",i, x[i], e[i]); } }

  • プログラムの処理について

    主クラスから別のクラスに送るときに 値一個ずつ6回送って処理した値を返すのと 値二個ずつ3回送って処理した値を返すのと どちらが早いのですか? 例1) main(){ int i,ans = 0; int num1[] = new int[6]; for(i=0;i<6;i++){ num1[i] = i+1; ans += ***.keisan(num[i]); } } class ***{ static int keisan(int num){ int x; x = num * 2; return x; } } 例2) main(){ int i,j,ans = 0; ***[] dt = new ***[3]; for(i=0,j=1;i<3;i++){ *** dt[i] = new ***(); dt[i].num1 = dt[i].num2 = 0; dt[i].num1 = j++; dt[i].num2 = j++; ans += ***.keisan(dt[i].num1,dt[i].num2); } } class ***{ public int num1; public int num2; static int keisan(int num1,int num2){ int x; x = (num1 + num2) * 2; return x; } } 上記の構文みたいな感じの状況です 構文が間違ってても気にしないでいただけるとありがたいです

    • ベストアンサー
    • Java
  • 数値解析の補間多項式

    (1)nを1以上の整数とし,X0,X1,,,Xnを相異なるn+1個の標本点とする。R上の関数f,g,hにおいて、gはfをX0,X1,,,Xn-1で補間し(つまり,g(Xi)=f(Xi),i=0,1,2,,,,n-1となる)、hはfをX1,,,Xnで補間するとき、関数    g(X)+(X0ーX)/(Xn-X0)×{g(X)ーh(X)} は、fをX0,X1,,,Xnで補間することを示したのですが質問があります。 まず補間するということはどんな意味を持っているのでしょうか?そしてこの問題の但し書きとしてf,g,hは多項式とは限らないとあったのですがではどう考えたらよいのでしょうか?? 最終的にどのように証明していけばよいかアドバイスお願いします★

  • 線形補間法プログラム(C++)

    C++言語で線形補間法のプログラムを組んで実行しているのですが、どしてもうまくいきません。ただ2倍の画像を作っているだけなのですが・・・。 以下プログラムを載せます、おおよその場所はわかるのですがどうすれば通るのかわかりません。どう直したらよいのか分かる方がいましたらご教授お願いします。 ※bmp[0]:現画像 pic:bmp[0]の縦横2倍の画像 // 線形補間法 // int zx = 2; int zy = 2; int i,j,m,n; float x,y,p,q; int xs = bmp[0]->Width/2; int ys = bmp[0]->Height/2; int d; for(i = -ys; i < ys; i++){ for(j = -xs; j < xs; j++){ y = i/zy; x = j/zx; if(y > 0){ m = (int)y; }else{ m = (int)(y-1);} if(x > 0){ n = (int)x; }else{ n = (int)(x-1);} q = y - m; p = x - n; if(q == 1){q = 0; m = m + 1;} if(p == 1){p = 0; n = n + 1;} if((m >= -ys)&&(m < ys)&&(n >= -xs)&&(n < xs)){ d = (int)((1.0 - q) * (1.0 - p) * (bmp[0]->GetPixel( m + ys, n + xs)) //おそらくこの辺に問題があるかと思われます。 + p * (bmp[0]->GetPixel( m +ys, n + xs)) + q * (1.0 - p) * (bmp[0]->GetPixel(m + 1 + ys, n + xs)) + p * (bmp[0]->GetPixel(m + 1 + ys, n + 1 + xs))); }else{ d = 0; } if(d < 0){d = 0;} if(d > 255){d = 255;} pic->SetPixel(i + ys, j + xs) = d; } } pictureBox2->Image = pic; }

  • このプログラムについて。

    #include <stdio.h> #include <stdlib.h> #include <math.h> #include <time.h> #include <string.h> #define NVALUE 30 #define MAXSIZE NVALUE struct City{ float x; float y; }; struct Country{ struct City cities[MAXSIZE]; int size; }; struct Itine{ float quality; int route[MAXSIZE]; int noc; }; float plength(struct Itine *tour,struct Country *pcountry) { int i,j; double dy; double dx; float length=0.0; for(i=0;i<tour->noc;i++){ j=i+1; if(j==tour->noc) j=0; struct City &c1=pcountry->cities[tour->route[i]-1]; struct City &c2=pcountry->cities[tour->route[j]-1]; double dx = c1.x-c2.x; double dy = c1.y-c2.y; length+=(float)sqrt(dx*dx+dy*dy); } return length; } なんか間違っていますか? エラーメッセージは、この部分 struct City &c1=pcountry->cities[tour->route[i]-1]; struct City &c2=pcountry->cities[tour->route[j]-1]; double dx = c1.x-c2.x; double dy = c1.y-c2.y; が、 ';'が'型'の前にありません と出ています。Visual Express 2005です。

  • このプログラムの問題点を教えてください。

    //数独 #include <iostream.h> const int num_Max = 100; //一辺のマス最高値 void input_data( int a, int [num_Max][num_Max] ); //入力関数 void yoko_data( int a, const int [100][100], int & ); //判定 void tate_data( int a, const int [100][100], int & ); void block_data( int a, int b, const int [100][100], int & ); main() { int num; //一辺のマスの数 int m; //一ブロックの一辺にあるマスの数 int okey; //numが正常か判別 int dx, dy, dz; int masu[num_Max][num_Max]; //全部のマス cout << "\n埋められた数独が正しいかどうか判断するプログラムです。\n"; while(1){ cout << "横何マスありますか? (4-100)>>" ; cin >> num ; if( num > num_Max || num < 4){ cout <<"範囲外です。再入力して下さい。" ; }else{ m = sqrt( num ); okey = num - m * m; if(okey != 0) cout << "その数字は数独として成り立ちません。再入力して下さい。"; } } //マスの入力 input_data( num, masu ); //横判定 yoko_data( num, masu, dx ); if(dx = 0){ //縦判定 tate_data( num, masu, dy ); if(dy = 0){ //マス判定 block_data( num, m, masu, dz ); } } if(dx == dy == dz == 0) cout <<"\n大正解♪ おめでっとー。\n"; return 0; } //入力 void input_data( int kazu, const int matrix[num_Max][num_Max]) { cout <<"\nマスの数字の入力を左上から順に、右へと行って下さい。"; for(int i = 0; i < kazu; i++){ for(int j = 0; j < kazu; j++){ cout << i+1 << "行目の" << j+1 << "列目 >>"; cin >> matrix[i][j] ; } } } //横列(行)を順に判定 void yoko_hantei( int kazu, int matrix[num_Max][num_Max]) { int kaburi; //かぶっているか判定する用 int vx = 0; //かぶっていたと判定されたかどうかを見る用 for(int i = 0; i < kazu; i++){ //横一列取り出しました。 for(int j = 0; j < kazu; j++){ for( int k = 1; j+k < kazu; k++){ kaburi = matrix[i][j] - matrix[i][j+k]; if(kaburi == 0){ //かぶってたらループから抜け出す vx++; break; } } } } if(vx > 0) cout << "\n横列(行)検索時に不適切な部分を発見しました。\n"; return vx; } //縦列(列)を順に判定 void tate_hantei( int kazu, int matrix[num_Max][num_Max]) { int kaburi; //かぶっているか判定する用 int vy = 0; //かぶっていたと判定されたかどうかを見る用 for(int j = 0; j < kazu; j++){ //縦一列取り出しました。 for(int i = 0; i < kazu; i++){ for( int k = 1; k < kazu-1; k++){ kaburi = matrix[i][j] - matrix[i+k][j]; if(kaburi == 0){ //かぶってたらループから抜け出す vy++; break; } } } } if(vy > 0) cout << "\n縦列(列)検索時に不適切な部分を発見しました。\n"; return vy; } //1ブロックごとに判定 void block_data( int kazu, int ruto, int matrix[num_Max][num_Max]) { int hantei[num_Max][num_Max] ; //判定する部分を切り出す用 int kaburi; int vz; //まずブロックごとに切り出してみる for(int i = 0; i < ruto-1; i++){ for(int j = 0; j < ruto-1; j++){ int h = 0; //何ブロック目か int g = 0; //そのブロックの何個目か (h-1)++; for(int l = 0; l < ruto-1; l++){ for(int k = 0; k < ruto-1; k++){ hantei[h][g++] = matrix[i * m + k][j * m + l]; } } } } //かぶっているか判定する for(int x = 0; x < kazu; x++){ for(int y = 0; y < kazu; y++){ for( int z = 1; z < kazu-1; z++){ kaburi = hantei[x][y] - matrix[x+z][y]; if(kaburi == 0){ //かぶってたらループから抜け出す vz++; break; } } } } return vz; } C++で作成したものです。 コンパイルエラーが出てしまうのですが、原因を教えていただけませんか?

専門家に質問してみよう