数学数学
正規化
次のサンプルは、データプロット(曲線)内のデータポイントを選択し、レイヤ内のすべての曲線をその点と同じ値に正規化します。 このコードは、複数曲線を持つグラフレイヤがアクティブであり、すべての曲線が同じX値を共有していることを前提としています。 この前提条件は、スペクトル分析ではよく行われることです。
GraphLayer gl = Project.ActiveLayer();
if( !gl )
return;
// Aある特定の曲線の1つの特定のポイントをクリックして選択
GetGraphPoints mypts;
mypts.SetFollowData(true);
mypts.GetPoints(1, gl);
vector vx, vy;
vector<int> vn;
if(mypts.GetData(vx, vy, vn) == 1)
{
// 選択したポイントのインデックスとy値を保存
int nxpicked = vn[0] - 1;
double dypicked = vy[0];
// レイヤ内のすべてのデータポイントをループ
foreach( DataPlot dp in gl.DataPlots )
{
// データ範囲、現在のプロットのy列
XYRange xy;
Column cy;
if(dp.GetDataRange(xy) && xy.GetYColumn(cy))
{
// y列からvector参照をy値に取得
vectorbase &vycurrent = cy.GetDataObject();
// vectorをスケールし、y値と選択したポイントを一致
vycurrent *= dypicked/vycurrent[nxpicked];
}
}
}
InterpolationExtrapolation補間/補外
ocmath_interpolate 関数は、線形、スプライン、Bスプラインのモードで補間/補外を実行するのに使われます。
// アクティブワークシートを4列にする
// 最初の2列に元のxyデータ
// 3列目に入力x値、4列目にy値
Worksheet wks = Project.ActiveLayer();
wks.SetSize(-1, 4);
DataRange drSource;
drSource.Add(wks, 0, "X"); // 1列目 - 元のxデータ
drSource.Add(wks, 1, "Y"); // 2列目 - 元のyデータ
vector vSrcx, vSrcy;
drSource.GetData(&vSrcx, 0);
drSource.GetData(&vSrcy, 1);
DataRange drOut;
drOut.Add(wks, 2, "X"); // 3列目 - 入力xデータ
drOut.Add(wks, 3, "Y"); // 4列目 - 補間したyデータ
vector vOutx, vOuty;
drOut.GetData(&vOutx, 0);
int nSrcSize = vSrcx.GetSize();
int nOutSize = vOutx.GetSize();
vOuty.SetSize(nOutSize);
int nMode = INTERP_TYPE_BSPLINE;
double dSmoothingFactor = 1;
int iRet = ocmath_interpolate(vOutx, vOuty, nOutSize, vSrcx, vSrcy, nSrcSize,
nMode, dSmoothingFactor);
drOut.SetData(&vOuty, &vOutx);
積分
Origin Cでは、積分を実行するのにNAGの積分ルーチンを使っています。 Origin CとNAGを使って、被積分関数、パラメータを持つ被積分関数、揺れのある被積分関数、不定積分、高次積分などで積分を実行できます。 次のサンプルは、NAGを使った積分を行う方法を示しています。
Origin Cコードには、NAG関数を呼ぶ前にNAGヘッダファイルをインクルードします。
#include <OC_nag.h> // NAGの宣言
簡単な積分関数
最初のサンプルは、1つの積分変数を持つ簡単な被積分関数で基本積分を実行する方法を示しています。
// NAG_CALL は適切な呼び出し表記。これを関数ポインタのように扱い
// 自分自身の被積分関数を定義
double NAG_CALL func(double x,Nag_User *comm)
{
int *use_comm = (int *)comm->p;
return (x*sin(x*30.0)/sqrt(1.0-x*x/(PI*PI*4.0)));
}
void nag_d01sjc_ex()
{
double a = 0.0;
double b = PI * 2.0; // 積分範囲
double epsabs, abserr, epsrel, result;
// 精度が十分でない場合、epsabs と epsrel を使って
// この値を目的の精度に拡張
epsabs = 0.0;
epsrel = 0.0001;
// 積分内の関数を評価するために部分間隔の最大数が必要。
// ほとんどの場合、200から500くらいが適切であり、お薦め。
int max_num_subint = 200;
Nag_QuadProgress qp;
NagError fail;
Nag_User comm;
static int use_comm[1] = {1};
comm.p = (Pointer)&use_comm;
d01sjc(func, a, b, epsabs, epsrel, max_num_subint, &result, &abserr,
&qp, &comm, &fail);
// 次の3つのエラー以外のエラーは
// 入力が不正か、割り当ての失敗。再び積分ルーチンを呼び出す前に
// メモリリークを避け、メモリ割り当てを
// 解放する必要がある
if (fail.code != NE_INT_ARG_LT && fail.code != NE_BAD_PARAM &&
fail.code != NE_ALLOC_FAIL)
{
NAG_FREE(qp.sub_int_beg_pts);
NAG_FREE(qp.sub_int_end_pts);
NAG_FREE(qp.sub_int_result);
NAG_FREE(qp.sub_int_error);
}
printf("%g\n", result);
}
パラメータ付きの積分関数
Integral, with Multiple Parameters
次のサンプルは、パラメータ付きの被積分関数で積分を定義し、実行する方法を示しています。 ユーザ定義の構造でパラメータが積分子に渡されます。 これは、スタティックな変数を被積分関数のパラメータとして使用し、スレッドを安全にします。
このサンプルは、NAGの不定積分子を使用することもできます。例えば、不定積分d01smc 関数を呼び出す行を有効にして、サンプルが不定積分を実行するのに使用します。
struct user // 被積分パラメータ
{
double A;
double Xc;
double W;
};
// ユーザによって供給される関数は、与えられたxで被積分関数の戻り値を返す
static double NAG_CALL f_callback(double x, Nag_User *comm)
{
struct user *param = (struct user *)(comm->p);
return param->A * exp(-2 * (x - param->Xc) * (x - param->Xc)
/ param->W / param->W) / (param->W * sqrt(PI / 2));
}
関数に対するパラメータをセットし、積分を実行する必要のある追加のパラメータを定義します。 そして、積分は、引数としてパラメータを渡し、1つの関数呼び出しで実行されます。
void nag_d01sjc_ex()
{
double a = 0.0;
double b = 2.0; // 積分間隔
// 次の変数は積分の精度を
// 制御するのに使用
double epsabs = 0.0; // 絶対精度、相対精度を使うには負値をセット
double epsrel = 0.0001; // 相対精度、絶対精度を使うには負値をセット
int max_num_subint = 200; // 最大部分間隔をセット、200-500がお勧め
// 結果は、アルゴリズムで返される近似積分値を保持
// abserrは |I - result| に対する上側境界の見積もり値
// ここで I は積分値
double result, abserr;
// Nag_QuadProgressの構造、
// max_num_subint要素を持つ内部的に割り当てられたポインタを含む
Nag_QuadProgress qp;
// NAGエラーパラメータ(構造)
NagError fail;
// パラメータはNAGの通信構造体で被積分関数に渡される
struct user param;
param.A = 1.0;
param.Xc = 0.0;
param.W = 1.0;
Nag_User comm;
comm.p = (Pointer)¶m;
// 積分実行
// NAGの不定積分子で使用可能な3種類の不定境界タイプ
// integrator Nag_LowerSemiInfinite, Nag_UpperSemiInfinite, Nag_Infinite
/*
d01smc(f_callback, Nag_LowerSemiInfinite, b, epsabs, epsrel, max_num_subint,
&result, &abserr, &qp, &comm, &fail);
*/
d01sjc(f_callback, a, b, epsabs, epsrel, max_num_subint,
&result, &abserr, &qp, &comm, &fail);
// エラーメッセージを出力してエラーをチェック
if (fail.code != NE_NOERROR)
printf("%s\n", fail.message);
// 次の3つのエラー以外のエラーは
// 入力が不正か、割り当ての失敗
// 再び積分ルーチンを呼び出す前にメモリリークを避け
// メモリ割り当てを解放する必要がある
if (fail.code != NE_INT_ARG_LT && fail.code != NE_BAD_PARAM
&& fail.code != NE_ALLOC_FAIL)
{
NAG_FREE(qp.sub_int_beg_pts);
NAG_FREE(qp.sub_int_end_pts);
NAG_FREE(qp.sub_int_result);
NAG_FREE(qp.sub_int_error);
}
printf("%g\n", result);
}
高次積分関数
2次以上の高次積分に対して、NAG積分子関数 d01wccを呼び出し、積分を実行します。
ユーザ定義のコールバック関数は、NAGのd01wcc関数に渡されます。
double NAG_CALL f_callback(int n, double* z, Nag_User *comm)
{
double tmp_pwr;
tmp_pwr = z[1]+1.0+z[3];
return z[0]*4.0*z[2]*z[2]*exp(z[0]*2.0*z[2])/(tmp_pwr*tmp_pwr);
}
メイン関数：
void nag_d01wcc_ex()
{
// 入力変数
int ndim = NDIM; // 積分の次元
double a[4], b[4];
for(int ii=0; ii < 4; ++ii) // 積分範囲
{
a[ii] = 0.0;
b[ii] = 1.0;
}
int minpts = 0;
int maxpts = MAXPTS; // 関数評価の最大数
double eps = 0.0001; // 精度をセット
// 変数出力
double finval, acc;
Nag_User comm;
NagError fail;
d01wcc(ndim, f_callback, a, b, &minpts, maxpts, eps, &finval, &acc,
&comm, &fail);
if (fail.code != NE_NOERROR)
printf("%s\n", fail.message);
if (fail.code == NE_NOERROR || fail.code == NE_QUAD_MAX_INTEGRAND_EVAL)
{
printf("Requested accuracy =%12.2e\n", eps);
printf("Estimated value =%12.4f\n", finval);
printf("Estimated accuracy =%12.2e\n", acc);
}
}
微分
ocmath_derivative関数は、スムージング無しで単純な微分を計算するのに使用します。 関数は、以下に示すように、ocmath.hで宣言されます。
int ocmath_derivative(
const double* pXData, double* pYData, uint nSize, DWORD dwCntrl = 0);
関数は、すべての欠損値を無視し、データ点とその隣のデータ点間の2つの勾配の平均を取ることで微分を計算します。 dwCntrl引数が0というデフォルト値を使う場合、関数は、データの方向が変わるときに入力されます。
if( OE_NOERROR == ocmath_derivative(vx, vy, vx.GetSize()) )
out_str("successfully");
dwCntrlがDERV_PEAK_AS_ZEROにセットされると、データの方向が変わった場合に関数に0が入力されます。
if( OE_NOERROR == ocmath_derivative(vx, vy, vx.GetSize(), DERV_PEAK_AS_ZERO) )
out_str("successfully");
|