Ajuste de datos a funciones trigonométricasHola,
Estoy luchando para ajustar los datos medidos a una función trigonométrica Ksin (ax + b) + c. esta es la naturaleza de mi proceso y estoy buscando una forma automática de ajustarlo mediante el código jsl.
Escribí un código que normaliza los datos para que la oscilación esté entre -1 y 1.
Al final, me gustaría obtener las constantes ayb (K y c también, pero es mucho más fácil)
ejemplo adjunto para mi tipo de datos.
BR,
Yoav
AnalíticaTabla de datosOtroScriptingRe: Adaptación de datos a funciones trigonométricasNo sé si esto es lo que quieres, pero puedes resolverlo usando la plataforma no lineal.Se adjuntan sus datos, junto con una columna de predicción (usando los parámetros para estimar) y un script que resuelve esos parámetros.Debo admitir que me tomó bastante tiempo obtener las suposiciones iniciales de los parámetros que permitieron que la solución convergiera (probablemente debido a mi falta de familiaridad con la función sinusoidal; ha pasado mucho tiempo desde que trabajé con tales funciones).Pero, conceptualmente, la plataforma no lineal logrará lo que desea.
Re: Adaptación de datos a funciones trigonométricasA continuación, se explica cómo obtener una estimación del a parámetro utilizando la función FFT. Los comentarios sobre segundos y horas suponen que los datos se recopilaron aproximadamente cada 4 segundos durante aproximadamente una hora; a la FFT no le importan las unidades. La FFT agrupa las frecuencias en los datos de manera que el primer intervalo después del intervalo de compensación de CC representa la frecuencia más baja que puede hacer un ciclo completo en los datos de entrada, el siguiente intervalo es de 2 ciclos completos, ... el intervalo del medio es ndata / 2 ciclos. Los contenedores después de que son un espejo del primero.
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);
Una respuesta fuerte en el cuarto contenedor después del contenedor de CC
Componente de frecuencia más grande a 0.00116448326055313 Hz o 859 segundos por ciclo
a_from_fft = 0.007316664113164;
El componente DC está relacionado con la media de los datos, la C valor.
Re: Adaptación de datos a funciones trigonométricasCreo que ayuda usar lo que llamo la forma estándar de una función seno: Y = K * sin (A * (x -
) + C.
En tus datos:
K = 1/2 del rango de valores
A = 2 * pi / período
B = desplazamiento a la derecha desde 0 hasta el punto de inicio de la onda sinusoidal
C = media de sus valores de Y
Las estimaciones iniciales para la regresión no lineal se pueden encontrar utilizando estas definiciones de la plataforma de distribución y el generador de gráficos:
K = 1/2 * (Máx. - Mín.) = 1/2 * (47 - 17,6) = 14,7
A = 2 * 3,14 / 880 = 0,007
B = 662
C = 34,1


El nonlinea r retroceder ion converge muy rápidamente usando este enfoque.El modelo convergió en 5 iteraciones utilizando estos valores. Vea abajo

Re: Adaptación de datos a funciones trigonométricasGenial, esa es exactamente la intuición que me faltaba.En cambio, manipulé los valores de los parámetros en la plataforma no lineal hasta que el gráfico se veía cerca de los datos, pero saber qué representan los parámetros (y obtener buenas suposiciones iniciales de los datos) es una forma mucho mejor de hacerlo.
Re: Adaptación de datos a funciones trigonométricasHola, gracias a todos por el tiempo que dedican. Este parece un gran enfoque para mi aplicación. pero, ¿cómo obtuviste "automáticamente" la duración del período? No entendí eso ...
Yoav
Re: Adaptación de datos a funciones trigonométricasPodrías usar un detector de cruce por cero , que sería un poco de JSL para revisar sus datos y notar cada vez que los datos pasan de negativos a positivos. Ese es el comienzo de un ciclo de onda sinusoidal. La distancia hasta el siguiente cruce por cero es el período.
Es posible que tenga problemas con datos ruidosos que parecen cruzarse varias veces.
El cruce por cero es más fácil de detectar que los picos.
Usé esa idea en Falla .