こんにちは、
測定データを三角関数のKsin(ax + b)+ c関数に適合させるのに苦労しています。これが私のプロセスの性質であり、jslコードでそれを自動的に適合させる方法を探しています。
振動が-1から1の間になるように、データを正規化するコードを作成しました。
最後に、aとbの定数を取得したいと思います(Kとcも同様ですが、はるかに簡単です)
私のデータ型の添付例。
BR、
ヨアヴ
これがあなたの望むものかどうかはわかりませんが、非線形プラットフォームを使用してこれを解決することができます。添付されているのは、予測列(パラメーターを使用して推定)およびそれらのパラメーターを解決するスクリプトとともに、データです。解が収束するのを可能にするパラメーターの最初の推測を得るのにかなりの時間がかかったことを認めます(おそらく正弦関数に精通していないためです-そのような関数を使ってから長い時間がかかりました)。しかし、概念的には、非線形プラットフォームはあなたが望むことを達成します。
見積もりを取得する方法は次のとおりですNS FFT関数を使用したパラメータ。 秒と時間に関するコメントは、データが約4秒ごとに約1時間収集されたと想定しています。 FFTはユニットを気にしません。 FFTは、DCオフセットビンの後の最初のビンが入力データで1つの完全なサイクルを作成できる最低周波数を表し、次のビンが2つの完全なサイクルであるように、データの周波数をビン化します...中央のビンはndata /です。 2サイクル。 後のビンは最初のビンのミラーです。
dt = Open( "z:\example data.jmp" );
// by hand, I made a formula column to look at x-lag(x) and
// decided the first 79 and last 8 were not representative...
rowstart = 80;
rowstop = 937;
timeStart = dt:x[rowstart];
timeStop = dt:x[rowstop];
ydata = dt[rowstart :: rowstop, {y}]; // throw away a few with irregular time
{r, i} = FFT( {ydata} );
nBins = N Items( r );
sampleRate = (timeStop - timeStart) / (nBins - 1); // imagine there were only 2 bins...thus -1.
// the samplerate is ~4 seconds/sample...
// the highest frequency this can resolve is 1 cycle/ 8 seconds. (.125 Hz)
// frequencies higher than that will cause aliasing unless they are filtered out before this.
// the bin [nBins/2] represents the maximum frequency that can be resolved
lastBinFrequency = 1 / (2 * sampleRate); // ~ 0.125 Hz
// the first bin [1] is the DC, or offset from 0, component. We'll ignore it.
// the next bin [2] represents the frequency that has one cycle over the length of the data.
// the lowest non-DC frequency is one cycle per hour (the length of time the data was collected).
binWidth = 1 / (timeStop - timeStart); // ~0.0003Hz, or ~1 cycle per hour
// note: lastBinFrequency == (nbins-1)*binWidth/2 -- each bin is this wide, in frequency
// the remaining bins are a mirror of the first and ignored.
// combine real and imaginary bin components to get an amplitude
p = Sqrt( r ^ 2 + i ^ 2 );
dtGraph2 = As Table( ((rowstart :: rowstop) - rowstart)` || p );
dtGraph2 << Graph Builder(
Size( 1300, 500 ),
Show Control Panel( 0 ),
Show Legend( 0 ),
Variables( X( :Col1 ), Y( :Col2 ) ),
Elements( Points( X, Y, Legend( 5 ) ) ),
SendToReport(
Dispatch( {}, "Col1", ScaleBox, {Format( "Best", 12 ), Min( 0 ), Max( nBins / 2 ), Inc( 100 ), Minor Ticks( 4 )} ),
Dispatch( {}, "Col2", ScaleBox, {Format( "Best", 12 ), Min( 0 ), Max( 30000 ), Inc( 5000 ), Minor Ticks( 1 )} ),
Dispatch( {}, "graph title", TextEditBox, {Set Text( "FFT Bins" )} ),
Dispatch( {}, "X title", TextEditBox, {Set Text( "bin number" )} ),
Dispatch( {}, "Y title", TextEditBox, {Set Text( "size" )} )
)
);
// there is a stong DC component in this data (all values > 0) so bin 1 is the max.
If( Loc Max( p ) != 1,
Throw( "something odd" )
);
p[1] = 0; // clear the DC value
// next max after DC:
maxBinIndex = Loc Max( p ); // 5. the data repeats about (5-1)==4 times in this frame.
Write(
"\!nBiggest frequency component at ",
binWidth * (maxBinIndex - 1),
" Hz or ",
Round( 1 / (binWidth * (maxBinIndex - 1)) ),
" seconds per cycle"
);
// where does a==.007 come from?
secondsPerCycle = 1 / (binWidth * (maxBinIndex - 1));
totalSeconds = timeStop-timeStart;
cyclesInSample = totalSeconds/secondsPerCycle;//about 4
timerange = timeStop-timeStart;
a_from_fft = cyclesInSample*2*pi()/timerange; // .0073...
show(a_from_fft);
DCビンの後の4番目のビンで1つの強い応答
0.00116448326055313Hzまたはサイクルあたり859秒での最大周波数成分
a_from_fft = 0.007316664113164;
DC成分は、データの平均に関連しています。 NS価値。
私が正弦関数の標準形式と呼んでいるものを使用すると役立つと思います:Y = K * sin(A *(x --B))+ C。
あなたのデータでは:
K =値の範囲の1/2
A = 2 * pi /期間
B = 0から正弦波の開始点まで右にシフト
C = Y値の平均
非線形回帰の開始推定値は、DistributionプラットフォームとGraphBuilderの次の定義を使用して見つけることができます。
K = 1/2 *(最大-最小)= 1/2 *(47-17.6)= 14.7
A = 2 * 3.14 / 880 = 0.007
B = 662
C = 34.1


ノンリネアr回帰このアプローチを使用すると、イオンは非常に迅速に収束します。モデルは、これらの値を使用して5回の反復で収束しました。 下記参照

素晴らしい、それはまさに私が欠けていた直感です。代わりに、グラフがデータに近くなるまで非線形プラットフォームでパラメーター値を操作しましたが、パラメーターが何を表しているかを知ること(およびデータから適切な初期推定を取得すること)は、はるかに優れた方法です。
こんにちは、あなたがそれに費やしてくれてありがとう。 これは私のアプリケーションにとって素晴らしいアプローチのように見えます。しかし、どのようにして期間の長さを「自動的に」取得したのですか? 聞き取れませんでした。
ヨアヴ
あなたは使用することができますゼロ交差検出器、これは、データをウォークスルーし、データが負から正に変わるたびに気付くためのJSLのほんの少しです。 これが正弦波サイクルの始まりです。 次のゼロ交差までの距離が期間です。
数回交差するように見えるノイズの多いデータに問題がある可能性があります。
ゼロ交差は、ピークよりも検出が容易です。
私はそのアイデアをグリッチ。