Matematica Interpolacion sen(x) / x o similares (osciloscopio)

Bueno, les cuento que pude arreglar un poco la etapa de entrada analogica y si, diversas resistencias que habia puesto en simulacion andaban bien pero en la realidad me hacian de pasabajos. Lastima que ahora que arregle esto, se me jodio el trigger, anda bien en determinadas condiciones. Mañana tocara revisar esa parte y luego a seguir con las matematicas!

Aca una cuadrada de 1mhz 5.44vpp medida con un osciloscopio hantek, tomada desde la entrada de mi osciloscopio (para contar la punta del mismo que esta distorsionando un poco por alguna razon)

1mhz-original.png

Y aca, mi osciloscopio midiendo la misma señal

1mhz-mio.png

Ahi la escala es 1us/div y estoy sampleandola a 20mhz (en vez de 40mhz) e interpolando 5 puntos en forma lineal. No adjunto la interpolacion con calculo splines debido a que empece a tener un error de desborde de vectores que tengo que corregir. De todos modos por lo que pude ver splines anda muy bien para cuadradas tambien.

Notese que en la parte de arriba esta bien plano, debido a que la señal se esta pasando de la escala vertical max y actuan los diodos de recorte. Pero noten la parte inferior como es muy parecida al osciloscopio comercial.

Tengo cada vez mas ganas de terminar este bicho de una vez para embarcarme a hacer uno portatil, mucho mas sencillo, que solo use el ADC del pic que es bien lento y un buffercito de entrada. Con la implementacion de splines de numeros enteros que desarrollo palurdo se podria obtener algo bastante interesante para ser chiquito, simple y portatil (con un GLCD como tantos osciloscopios que andan dando vueltas por ahi)

Ah! Adjunto tambien una simulacion en proteus de la entrada (1 solo canal + trigger), les pido disculpas por lo horriblemente chancha que esta su distribucion. Aclaro que no estoy usando en la realidad el mismo operacional por lo tanto, en el aparato, el preset de entrada esta hacia VEE en vez de VCC y el valor no es el mismo sino mas grande.
 

Adjuntos

  • entrada_mux_selector-rail-to-rail.zip
    23.2 KB · Visitas: 7
Tengo cada vez mas ganas de terminar este bicho de una vez para embarcarme a hacer uno portatil, mucho mas sencillo, que solo use el ADC del pic que es bien lento y un buffercito de entrada. Con la implementacion de splines de numeros enteros que desarrollo palurdo se podria obtener algo bastante interesante para ser chiquito, simple y portatil (con un GLCD como tantos osciloscopios que andan dando vueltas por ahi)

Ojo que los GLCD son en general bastante lentos.
He trabajado con algunos de 128x128 y 128x64 pixels y tienen cierto tiempo de refresco sólo
tolerable por los que alguna vez sufrimos las XT 286, como yo. :D

Hace un tiempo estoy poniendole mis horas sueltas a un aparatito construido con tu misma filosofia.
No había pensado trabajarle el interpolado porque no hay mucho lugar dónde poner más puntos.
Sacando algunas columnas para textos me quedan unos escasos 100 pixels de ancho para la gráfica.
A ver si les gusta la pinta.
La placa de adaptación es externa para poder probar cosas sin tener que desarmar todo.
En este caso toma señal analógica de dos "electret".
 

Adjuntos

  • ODM2B.JPG
    ODM2B.JPG
    46.7 KB · Visitas: 28
  • ODM2B_2.JPG
    ODM2B_2.JPG
    28.8 KB · Visitas: 23
Última edición:
Yo empece con una TI99 asi que estoy acostumbrado a la lentitud jaja. Otra idea un poco mas complicada pero aun factible es usar una pantalla de notebook vieja. Venden por ebay una plaquita que es una interfaz vga-pantalla lcd con el elevador de tension y todo. Si con un micro extra (tipo controlador de pantalla) se genera la señal vga se podria hacer un osciloscopio digital de mesa. Las posibilidades estan ahi pero primero voy a terminar este que ya esta casi para ir a un gabinete.
 
las pantallas de 7" de las tablets casi son mas baratas que los glcd. Otra cosa es con que controlarlas... Hay que pensarlo muy bien a ver si va a salir mas caro el collar que el perro jejejeje
 
Código:
function [xmax ymax xmin ymin] = maxminsplines(samples)

  Y_1=samples(length(samples));

xmax=0;
ymax=0;
xmin=0;
ymin=0;

  for i=1:(length(samples)-1)

	Y0=samples(i);
	Y1=samples(i+1);
	if (i==length(samples)-1) 
 	  Y2=samples(1);
	else 
        	  Y2=samples(i+2);
	end

	S=Y2+Y0-2*Y1;
	T=Y1+Y_1-2*Y0;

	A=(S-T)/2;
	B=(2*T-S)/2;
	C=(Y1-Y_1)/2;
	D=Y0;
                  
        if  A==0
		X1=-C/(2*B);
		X2=X1;
	else
		X1=(-2*B+sqrt(4*B^2-12*A*C))/(6*A);
		X2=(-2*B-sqrt(4*B^2-12*A*C))/(6*A);
	end

	Y1=A*X1^3+B*X1^2+C*X1+D;
	Y2=A*X2^3+B*X2^2+C*X2+D;

	if X1>0 && X1<1 && imag(X1)==0
		if Y1 >ymax
			ymax=Y1;
			xmax=i+X1;
		end
		if Y1 <ymin
			ymin=Y1;
			xmin=i+X1;
		end
	end

	if X2>0 && X2<1 && imag(X2)==0
		if Y2 >ymax
			ymax=Y2;
			xmax=i+X2;
		end
		if Y2 <ymin
			ymin=Y2;
			xmin=i+X2;
		end
	end

	if D>ymax
		ymax=D;
		xmax=i;
	end 
	if (A+B+C+D)>ymax
		ymax=A+B+C+D;
		xmax=i+1;
	end
	if D<ymin
		ymin=D;
		xmin=i;
	end 
	if (A+B+C+D)<ymin
		ymin=A+B+C+D;
		xmin=i+1;
	end 
	Y_1=Y0;
	end
 return
end

Palurdo, he implementado en mi programa la funcion para extraer datos de los samples crudos que me detallas mas arriba. Los resultados parecen ir teniendo forma.

Tuve que establecer el valor inicial de xmin = 9999 ymin = 256 (mi maximo es 255) para que esos datos tomen algun valor sino nunca entra a los IF. Por lo que veo en tu codigo los inicializas a cero, entonces supongo que tus "samples" incluyen numeros negativos (en mi caso, de 0-255)

Los valores que obtengo para YMax YMin se ven muy razonables y varian cuando cambio la amplitud de la señal, me falta traducirlos a volts segun la escala que este la entrada.

Sin embargo, para Xmax y Xmin obtengo datos un poco erraticos, esta fue mi implementacion para calcular frecuencia (pongo unidades aqui en la formula asi se entiende que es cada parte)

sampling_rate = [1=25ns,2=50ns,3=100ns]

dDiff = ((Abs(xmax - xmin) * 2 [dobleperiodo] ) * (25ns * sampling_rate)) / 1000000000 [cuantos ns es 1 segundo]

dDiff = 1 / dDiff

La frecuencia de la señal de prueba es otra vez 20khz y obtengo, cuando 1 solo periodo de la señal ocupa toda la pantalla (todos los samples) valores similares a 20.000hz pero que varian bastante rapido entre 19 y algo y 21 y algo...

Sin embargo, en escalas mayores donde tengo varios picos de la onda seno x pantalla, la frecuencia calculada (por ende, los valores xmax y xmin) pegan saltos muy grandes y muy rapidos.

Una aclaracion importante seria que en tu funcion existe esta linea:

if X1>0 && X1<1 && imag(X1)==0

Que yo la tuve que cambiar (para probar nomas) por esto:

if X1>0 && X1<1

Debido a que no se que hace la funcion imag(x) (sera numeros imaginarios o algo asi?)

Si pudieses aclararme que hace, ya que no esta definida en el codigo y supongo que es alguna funcion de octave/mathlab, quiza solo eso sea el problema que estoy teniendo, ya que sin esa ultima condicion, en casi toda vuelta del ciclo for entro al IF.

Muchas gracias!
 
En Matlab imag(x) devuelve la parte imaginaria del número x.
Lo necesita para asegurarse que la x (que sale de la raíz de una ecuación cuadrática) sea un número real.

Lo raro es que en varios lugares usa: xmin = i + X2, y similares.
Los complejos no están ordenados, entonces no tiene sentido decir que un nro. complejo sea mínimo o máximo de un conjunto.
 
Última edición:
En Matlab imag(x) devuelve la parte imaginaria del número x.
Lo necesita para asegurarse que la x (que sale de la raíz de una ecuación cuadrática) sea un número real.

Lo raro es que en varios lugares usa: xmin = i + X2, y similares.
Los complejos no están ordenados, entonces no tiene sentido decir que un nro. complejo sea mínimo o máximo de un conjunto.

Mea culpa. Ya sabes que en ingenieria usamos j en lugar de i para variable compleja asi que no meimporto usar i como varianble de indice. Igual habiendo puesto xmin=ind+X2 se habria visto mas claro. De todas formas me aseguro que todos los valores sean reales.

Seaarg:

En tu caso, simplemente asegurate de que lo que hay dentro de la raiz cuadrada sea mayor que 0. En caso de queno sea asi no calcules X1, X2, Y1, Y2 porque esos valores no caen dentro de el eje real.



la funcion con inicializar xmin, xmax ymin e ymax con x e y del primer sample debe ser suficiente.
 
Última edición:
Bueno, en vb.net tenemos un constructor de clase numerica tipo "complex" al cual asigne (luego de lo que me explico asherar) el resultado de la operacion (defini variables X1 y X2 como complejos) y en todas los calculos donde intervienen esas variables, uso la propiedad real, que devuelve la parte real de un numero complejo. En el IF comparo la propiedad "imaginary" contra cero. Creo que logre hacer lo mismo que tenes hecho en tu codigo con estas herramientas.

Sin embargo el problema persiste y creo que encontre la causa aunque no sabria como solucionarlo:

Probe implementar tanto tu funcion, asi como tambien una funcion mia bien "bruta" donde simplemente promedio 3 valores continuos de Y en el tiempo, buscando maximo y minimo a lo largo de todos los samples. Y el resultado de ambas es similar en comportamiento (no asi en valores, estoy 99% seguro que los valores devueltos por tu funcion son mas acertados)

Resulta que, si tengo varios picos de la onda seno en pantalla (o sea, en el set de samples), la funcion va recorriendolos, buscando el maximo y minimo valor de Y (luego de los calculos claro) y sucede lo siguiente:

Para un set de 896 datos: 0-895 en eje x

- En una primera pagina, digamos que encuentra el minimo Y en el sample 80 y el maximo en el 160
- En una segunda pagina, encuentra el minimo en el mismo sample 80 pero el maximo por el 500 (debido a que el "maximo" del punto 160 no fue mayor a un pico que hubo en esta pagina, en el sample nro 500)
- En una tercera pagina, el minimo lo encuentra en 200 y el maximo en 896

Por supuesto los valores que pongo son de ejemplo extremo de lo que sucede. La cuestion es que claro, como la señal no es algo perfecto y una cresta puede tener 1 punto mas que otra tanto para arriba como para abajo, detecta el periodo entre 1 cresta y otra, o entre 2 crestas, o 3 y asi, haciendo que el calculo de distancia entre crestas que hago sea variante entre cada "pagina" de samples.

Intente hacer una especie de atajo para evitar esto, reemplazando:

Código:
if Y1 >ymax 
  ymax=Y1;
  xmax=i+X1;
end

por:

Código:
if Y1 >ymax + 10
  ymax=Y1;
  xmax=i+X1;
end

para que asi, solo cuando el valor es bastante mas grande cambie la variable ymax (seria como una deteccion de pico maximo y minimo digamos) con iguales resultados.

Tambien probe implementar un detector para que solo cuando la señal "sube" compare YActual con Ymax y solo cuando "baja" Yactual con Ymin con los mismos resultados (claro! me di cuenta despues de probar que es lo mismo)

Los valores de YMAX e YMIN se mantienen mas o menos constantes, con variaciones de 1 o 2 enteros nomas, por lo que el calculo de vmin y vmax de la señal ya estoy en condiciones de hacerlo pero el calculo de la frecuencia se me resiste por esa variabilidad.

Estoy pensando como puedo hacer para detectar 1 cresta y su valle siguiente, o viceversa, para capturar solamente 1 periodo (siempre duplicando la distancia en X claro porque la F es entre cresta y cresta)

Esto se pone cada vez mas lindo! :)

Edito con una pregunta:

¿Que opinan de esto? Si en la funcion de palurdo pongo mas codigo que haga algo como:

- Detectar todas las crestas, es decir, para cada vez que y > ymax almaceno ese valor como ymax pero, una vez que el valor de Y sea menor a ymax (con cierto margen, o sea, cuando empieza a bajar) "reseteo" ymax para detectar la siguiente cresta... y asi sucesivamente hasta terminar los samples.

Luego, teniendo las varias distancias en X entre las distintas crestas, hacer un promedio entre todas ellas y que eso sea la diferencia xmax - xmin.

Igualmente, la deteccion de que termino una cresta es mas facil decirla que programarla.
 
Última edición:
Claro, es que la función detecta el máximo y mínimo global y lo que necesitas es detectar 2 máximos locales o dos mínimos locales. Se puede modificar la función para que además de los globales te detecte todos los locales de la señal. Como no sabemos cuantos máximo y mínimos locales puede tener la señal en la gráfica, tendrás que ir guardándolos en un vector de máximos y un vector de mínimos.

Bueno, pues cada vez que se entra en uno de los if X1>0 y X1<1 entonces ese X1 (o X2) es o bien un máximo o un mínimo. Para saber si es un máximo o un mínimo haces el cálculo de 6*A*X1+2*B y si ese resultado es negativo tienes un máximo, y si es positivo tienes un minimo, y si es 0 no es ni un máximo ni un mínimo sino un punto de inflexión (que a ti no te es útil). Además X es un máximo local si X[i-1] y X[i+1] son menores, y es un mínimo local si son mayores.

Por lo tanto por cada dos puntos consecutivos yo comprobaría si cada uno de los puntos X, X1 y X2 son máximos o mínimos y en caso de que algunos lo sean, añadirlos a un vector.

Por definición, los puntos inicio y final no pueden ser máximos locales, por lo que no te fijes si X[0] y X[895] con mínimos y máximos locales.


Por otro lado este método sólo sirve para calcular la frecuencia de una senoidal pura o de una onda monótona (cuyas crestas y valles coinciden con su periodo), pero falla cuando la señal sea la suma de varias senoidales, por ejemplo, porque calcularás la frecuencia de un armónico de la señal. Lo que puedes hacer es calcular amplitud y frecuenci entre un máximo y un mínimo local, despues generar una senoidal con esa amplitud y frecuencia (y ajustando la fase para que coincida el máximo y el mínimo de la senoidal generada con los máximos y mínimos obtenidos) y le restas esa senoidal a tu onda original. repites el proceso, y si la amplitud de la señal resultante es menor que la de tu senoidal, pues ya tienes la componente principal. En caso contrario vuelves a calcular máximos y mínimos y vuelves a generar otra senoidal y la restas, y así hasta que la senoidal que generes sea mayor que la amplitud de la onda que te queda. En ese caso de la lista de senoidales generadas, la que tenga mayor amplitud es la componente principal de frecuencia.
 
... Otra idea un poco mas complicada pero aun factible es usar una pantalla de notebook vieja. Venden por ebay una plaquita que es una interfaz vga-pantalla lcd con el elevador de tension y todo. Si con un micro extra (tipo controlador de pantalla) se genera la señal vga se podria hacer un osciloscopio digital de mesa. .
Estuve gugleando y encontré formas simples de acoplar salidas VGA con entradas a LCD ( ENLACE al artículo 1) y también con RGB+Sync tipo salida de PlayStation (ENLACE al artículo 2).
Aparentemente todo pasa por tener el pinout del LCD. Si no ... (n)

laptop-lcd-display-to-vga-interface3.thumbnail.jpg


Yo tengo una Sharp LM130SS1T611, que venía en una Compaq Presario, hoy fenecida de muerte natural.
Si alguien tiene el pinout por favor le pido que lo suba por acá (y si no es demasiada molestia me avise por MP).
Mil gracias.
 
Última edición:
Buscando programas que implementen funciones e interpolaciones, he caido en esto:
http://jonw0224.tripod.com/ppmscope.html

Se lo ve prometedor al chiquitin. El soft es bastante completo y tiene interpolaciones splines y sinc

Viene con soft para PC, el firmware del pic (en assembler!) PCB, lista de componentes, etc.etc.

Editado: Pongo este post en el thread de asherar sobre osciloscopios, ya que me parece mas adecuando que en este. Perdon!
 
Última edición:
Les cuento, al final pude obtener la frecuencia de la señal con un metodo medio feo que es un detector de cruce por cero. Cuento los samples entre 2 cruces por cero, eso lo multiplico por 2 y con el valor de tiempo x sample obtengo el periodo y de ahi la frecuencia.

No es un buen metodo, si es simple.

Adapte desde codigo C hacia vb.net unos 3 o 4 codigos que pude encontrar para calcular FFT, y al final me quede con el mas prolijo y completo, que es el que adjunto (el archivo se llama func.c) y que pertenece al soft del osciloscopio PPMScope

Creo haberlo adaptado bien, reemplazando los tipos de datos, las copias de vectores, los operadores especificos de C, etc.

Para el que lo quiera ver, la cosa esta en la funcion sigFreq() donde se aplica una ventana de Hann a los samples, luego se realiza la FFT y luego se busca la frecuencia fundamental pasandole las magnitudes calculadas en la FFT a traves de fundFreqBin, que por lo que entendi, calcula la frecuencia dividiendo la FFT en contenedores que representan cada frecuencia posible.

Dicho esto: NO HAY FORMA! algo esta mal en mi codigo o en mi interpretacion de los resultados, porque SIEMPRE, sea la que sea la frecuencia que estoy sampleando, obtengo que la mayor magnitud esta en la posicion 256 de la FFT, que es la mitad de los datos enviados a la misma. Si uso la mitad de los datos de la FFT (n/2) entonces la posicion siempre es 128 y asi sucesivamente.

Lo mismo me paso con otros algoritmos para hacer FFT que fui traduciendo desde C a net.

¿Cual puede ser mi concepto errado por el cual esto no anda? Es decir: Puede que este interpretando mal la salida de la FFT en mi cabeza o habra un error de codigo?

Mi interpretacion es: La salida de la FFT es el dominio de las frecuencias. El valor de vector de mayor amplitud indica la posicion de la frecuencia fundamental. Luego: Como se cual es dicha frecuencia? (relacion con sampling rate 100% seguro pero no se como aplicarla) Esto esta definido en fundFreqBin() pero no usa sampling rate.

Aclaro que no pido que me resuelvan el programa, sino que me ayuden a entender el resultado de la FFT.

Adjunto tambien un archivo FFTAnalisis.zip que es un excel donde puse, en la columna A los samples originales (señal 20khz sampleados a 5.000.000 hz) y en la columna B un dump de lo que devuelve la funcion faFourier, en el vector de magnitudes.
 

Adjuntos

  • func.txt
    16.6 KB · Visitas: 10
  • FFTAnalisis.zip
    14.8 KB · Visitas: 15
Esta te devuelve el valor, escalado en frecuencia, del índice del máximo de la FFT.
/* Returns the fundamental frequency of the signal
* arr[]: array in which to find the fundamental frequency
* len: length of the array
*/
double sigFreq(double arr[], int len)
{
...
}
Y usa ésta otra para buscar el índice.
/* A Helper function. (Una función de ayuda)
* Find the fundamental frequency bin, assumes no peak merging
(Encuentra el índice de la frecuencia fundamental, supone que no hay picos mezclados)
*/
int fundFreqBin(double magn[], int len)
{
...
}
Esta te devuelve el índice de la muestra de mayor amplitud.
Si el error es acá, debe ser que te queda siempre una muestra de valor muy alto al final, y por eso te queda ahí al salir.

Otra cosa más:
En la planilla hay un dato, en el casillero C1: (4,09E+75) que no me parece normal.
Hay otros similares por ahí, y por eso la gráfica te da las abscisas entre 0 y 10^80.
La forma de la FFT tiene muchos picos de ese orden de magnitud.
Si es la Transformada de la sinusoide de la primera gráfica, debería ser un solo pico en la f del seno.

No probaste buscar algo ya hecho para PIC en la página de Microchip ?
Buscando aquí por ejemplo De paso ya nos queda para el portátil !!!

En el foro se ha tratado el tema: https://www.forosdeelectronica.com/f24/fft-pic-13856/

Implementation of Fast Fourier Transforms
pdf.gif


Reading and Using Fast Fourier Transforms (FFT)
pdf.gif
 
Última edición:
Bien, primero que nada voy a revisar si la traduccion de C a vb.net la hice bien (creo que si) de esta funcion

Código:
int fundFreqBin(double magn[], int len)
{
    double avg, a;
    int i;

    avg = sigMax(&magn[1], len-1)/2;

    i = 0;
    a = magn[i];
    while(magn[++i] < a)
        a = magn[i];
    a = magn[i];
    while(magn[++i] > a || a < avg)
        a = magn[i];
    i--;

    if(i >= len - 1)
         return 0;
    return i;
}

A esto en vb.net

Código:
    Private Function fundFreqBin(ByVal magn() As Double, ByVal len As Integer) As Integer

        Dim avg, a As Double
        Dim i As Integer
        Dim tmagn(magn.Length - 2) As Double
        Array.ConstrainedCopy(magn, 1, tmagn, 0, len - 1)
        avg = sigMax(tmagn, len - 1) / 2

        i = 0
        a = magn(i)
        i = 1   '++i del while
        While magn(i) < a   ' La primera vez = 1
            a = magn(i)
            i = i + 1
        End While
        a = magn(i)
        i = i + 1
        While magn(i) > a Or a < avg
            a = magn(i)
            i = i + 1
        End While
        i = i - 1
        If i >= len - 1 Then
            Return 0
        Else
            Return i
        End If

    End Function

Luego, tengo un problemon y es el siguiente:

BIN = Fs/N = 40000000/1024 = 39062Hz.

Eso significa que, para 1024 muestras en la maxima frecuencia de sampleo del osciloscopio, la variacion entre 1 bin y otro es ENORME!

Se me ocurre que podria tomar, por ejemplo 1 muestra de cada 2 para asi "bajar" la F de sample a la mitad y completar el resto con ceros o algo por el estilo pero aun asi los saltos entre BIN son grandes.

Tambien sino, fabricar muchas mas muestras usando la interpolacion splines de palurdo para achicar de esta forma el bin size. El problema de esta forma seria que tendria que procesar una FFT para, digamos, 1024*10 muestras, con todo lo que eso significa en consumo de memoria y en velocidad de proceso (ya hasta la PC se pondria lenta creo)

Es decir: para hacer FFT me encuentro con no solo los problemas tecnicos de programacion, etc, sino de que la Fs del osciloscopio es demasiado grande en la escala de tiempos mas chica. ¿Es correcta esta interpretacion?
 
A mí la traducción de las rutinas me parece correcta.

Luego, tengo un problemon y es el siguiente:

BIN = Fs/N = 40000000/1024 = 39062Hz.

Eso significa que, para 1024 muestras en la maxima frecuencia de sampleo del osciloscopio, la variacion entre 1 bin y otro es ENORME!

Se me ocurre que podria tomar, por ejemplo 1 muestra de cada 2 para asi "bajar" la F de sample a la mitad y completar el resto con ceros o algo por el estilo pero aun asi los saltos entre BIN son grandes.

Tambien sino, fabricar muchas mas muestras usando la interpolacion splines de palurdo para achicar de esta forma el bin size. El problema de esta forma seria que tendria que procesar una FFT para, digamos, 1024*10 muestras, con todo lo que eso significa en consumo de memoria y en velocidad de proceso (ya hasta la PC se pondria lenta creo)

Es decir: para hacer FFT me encuentro con no solo los problemas tecnicos de programacion, etc, sino de que la Fs del osciloscopio es demasiado grande en la escala de tiempos mas chica. ¿Es correcta esta interpretacion?

No entiendo bien tu problema.
Si no querés tener un BIN demasiado grande bastaría con expresar tus frecuencias en kHz o en MHz. La cantidad representa lo mismo pero con un número más bajo.
Si usas los splines en realidad te aumenta el valor en lugar de bajar.
Creo que por ese tema algunas rutinas hacen cuentas con el índice y luego calcula la frecuencia en las unidades que corresponde.
Pero no sé para qué usas luego esos números, así que tal vez mi observación no tenga sentido.
 
Última edición:
Tal vez estoy interpretando mal los calculos que se hacen con la FFT, tengo entendido lo siguiente:

- Tengo un vector de 1024 samples en dominio de tiempo
- Paso esos samples por FFT
- Tengo ese vector, pero en dominio de frecuencia
- Aplico calculo de magnitud para cada elemento del vector
- Obtengo el indice del vector donde se encuentre la mayor magnitud, salteando el x[0] que es DC
- Calculo el bin size con la Fs de esta forma: bin = Fs / 1024
- Multiplico el indice de mayor magnitud por este bin size para obtener la F fundamental de los samples.

Entonces, para Fs = 40000000 en esa cantidad de samples el bin es = 39062,5

Supongamos que la mayor magnitud de la FFT esta en el indice 1, entonces el resultado del calculo de frecuencia fundamental va a dar 39062,5 HZ (y de ahi para arriba, asi que algo esta mal!! no podria obtener una F menor a 40khz, para esta cantidad de samples, a este sample rate tan grande)

Llego a esa conclusion porque en diversas pruebas, con diversos set de samples, estoy obteniendo esos indicios, pero estoy casi seguro de que hay algo mal o en las presuposiciones, o en las rutinas que hice.

Consegui un programita en vb 6.0 que hace el grafico de la FFT de una senoidal generada por el mismo (con funcion sin()) y el mismo devuelve en que indice esta la max amplitud...

... entonces, le reemplace la generacion artificial de la senoidal por un conjunto de mis samples:

Resulta que este si me da el indice siempre en los primeros lugares (no como mi rutina que siempre me daba en 256 o similares) y, como sample esas pruebas a 20us/div, hice la siguiente cuenta:

bin = 5000000 / 1024 = 4882,81

La aplique y parece ir funcionando, me da valores en HZ similares a la senoidal que yo genero con la placa de sonido. Por supuesto, la resolucion es de 4.8khz entre un indice y el siguiente, por lo tanto, la resolucion es pesima pero al menos es logica.

Tal vez aqui tenga que releer algo que explico palurdo mucho mas arriba, sobre la frecuencia fundamental que cae entre dos indices enteros.

Es curioso tambien como, en este programita esta la opcion de usar 256, 512, 1024, 2048 y 4096 samples. En 1024 (yo tengo 896 muestras y completo con numeros 127 que seria el GND) anda bastante bien, pero a medida que aumento la frecuencia de la senoidal sampleada, tengo que subir la cantidad de samples. Sin embargo, si uso la opcion de 4096 samples para todas las pruebas, en la mayoria de las situaciones la magnitud maxima de la FFT se encuentra en el indice 0 (como si fuera DC)

En fin, me voy a volver loco jeje. Lo que voy a hacer ahora es estudiar otras implementaciones de la FFT, en php por ejemplo asi puedo ver si estoy cometiendo un error y donde.

Me encantaria poder interpretar los simbolos matematicos de la descripcion de la FFT. Solo conozco el de sumatoria :oops: entonces no puedo hacer un algoritmo en codigo directamente a partir de la formula. (ay, la ignorancia!!!)

Edito: Me encontre un error en mi traspaso de codigo desde C a vb.net y ahora en vez de tener la mayor amplitud en la mitad de la FFT siempre, la tengo siempre en el indice de la FFT = 1 (descartando el 0 que es DC)

Tambien adapte ese programita vb 6.0 que comente que funciona, y le puse que trabaje con 8192 samples (896 reales entre -127 y 127 y padding a 0). Desde frecuencias en el generador de algunos cientos de hz hasta los 20khz me da una medida que parece ser buena, para una Fs de 5.000.000

Obviamente, este programa que encontre una un algoritmo un poco distinto para hacer la FFT. Podria adaptarlo pero el codigo C que tenia me parecia mas elegante (este codigo de vb 6 funciona pero es bastante chancho a mi gusto)

Quisiera aclarar que, para el calculo de frecuencia deje de intentar aplicar max y min locales como me decia palurdo, y me emperre en sacar adelante FFT, por la razon de que luego, usando la misma FFT puedo sacar otros valores como THD, etc. Y aclarar que la funcion que me brindo palurdo quedo super firme para vmax y vmin. Esos estan saliendo perfectos.

Sigo intentando buscarle la vuelta. Les comentare resultados!
 
Última edición:
Vuelvo con mas novedades. Pude implementar la FFT :D

Consegui un programita hecho en vb.net que hace la FFT, tiene aplicacion de ventana de hamming, hanning y blackman (o ninguna) y funcionaba bastante bien, asi que agarre el codigo fuente y lo adapte a mi programa.

Estoy usando la barbaridad de 131072 muestras (896 de verdad, el resto cero) y es demasiado rapido, las procesa sin problemas!

Se eligio semejante numero para poder disponer de un bin size de 40000000 / 131072 = 305hz

Lo proximo seria agregarle un pequeño codigo que cambie la cantidad de samples que se le envian a la FFT segun el sample rate del ADC, porque es inutil mandar tantos cuando sampleo mas lento que 40mhz. De esa forma, seria mas rapido aun.

El algoritmo es bien feo en cuanto a codigo de programacion (uso de goto) :) pero lo importante es que funciona muy bien. Lo adjunto aqui por si alguien lo quiere investigar, yo solamente le tuve que adaptarlo a mi señal, sampling rate variable, etc y hacer el calculo de magnitud sqrt(re^2+im^2) para poder obtener la frecuencia fundamental.

Luego de unas cuantas adaptaciones, esto es lo que logre visualizar:

confft.png

Ahi se ve la señal de los 2 canales, la frecuencia y periodo (aun por adaptar a algo legible) de los mismos, los niveles de señal obtenidos por splines, la ventana de grafico de FFT para el canal A (verde) mostrando el espectro, y deje tambien ahi la ventanita del programa generador de señal, para comparar lo que dice que esta generando con lo que arrojan de resultado las matematicas.

En fin, parece que todo va bien! ahora ya me puse ambicioso y le quiero sacar mas informacion a la FFT, como THD, ruido, etc.

Si alguien quiere el programa completo del osciloscopio, digame y con gusto lo subo aqui. Con explicaciones del codigo mediante, otra persona puede hacerse su placa de captura y utilizar este soft (que aun le falta mucho! esta verde, pero pinta bien)
 

Adjuntos

  • Fft2.rar
    36 KB · Visitas: 11
Excelente proyecto, me están dando ganas de armarme uno, voy a buscar por eBay algunos ADCs rápidos.

PD: Con solo un "poco" mas de trabajo lo podrías vender por medio millón de dolares, solo necesitas ampliar el ancho de banda a 62GHz y tomar 160GS/s. Vídeo de alguien "jugando" con un osciloscopio de medio millón de dolares
 
y en que quedo el proyecto?? ya funciono???

Si que funciono! luego de aprender mucho, renegar mucho y recibir mucha ayuda, ahora tengo una herramienta bastante util para proximos diseños :D

Aca va una imagen de una cuadrada generada por un cristal de 8mhz y dividida con un contador. La señal es de 2mhz y esta siendo sampleada a 40mhz, utilizando la interpolacion SPLINES de palurdo (no me canso de agradecerte! jeje) Tambien se ve la ventanita de FFT donde veo la fundamental y las distintas armonicas. Dicha ventana la deje chica para que entre en esta captura de pantalla pero como tengo 2 monitores, normalmente arrastro esa ventana al 2do monitor y veo la FFT en full screen con todo su detalle.

2mhzcuadrada.png

En esta otra, vemos 2 senoidales en canal A y B, y ademas lo que se ve como "fantasma" por detras son 2 senoidales de la referencia. Le adicione al programa poder guardar una captura para despues cargarla de fondo como referencia visual (cuando comparo señales, por ejemplo, para ver si un filtro pasa bajos filtra mas o menos segun valores de componentes)

senoidal.png

Y por ultimo, pongo una foto del bicho en su gabinete temporal-definitivo jeje, son tapas de lectora de CD recicladas, ahora tengo que hacerle el frente nomas :LOL:

gabinete.jpg

Mejore varios aspectos del programa y ahora anda mucho mas rapido tambien!

Caracteristicas generales:

2 Canales + trigger externo
sample rate max: 40mhz
Escalas de tiempo: Desde 50ns/DIV a 1 S/DIV
Memoria: 896 bytes x2

Analogica:
Entrada de señal 1:1 = 10vpp (100vpp con el selector de la punta) ITEM A MEJORAR
Ancho de banda: 75-130mhz (no probado)
Ancho de banda probado: 4mhz (con señal cuadrada)
Divisores: 1:1, 2:1
Multiplicadores: 1x2, 1x10
AC/DC

Trigger:
Externo, Canal A, Canal B. Nivel seleccionable con mouse (5vpp a mejorar) Trigger automatico, single shot, slope + y -

Software:
Interpolacion SPLINES y Lineal para bases de tiempo 1us/DIV hacia abajo. Para arriba es captura real
Guarda y carga a memoria un DUMP de la señal
Guarda y carga a memoria una referencia de señal para mostrar de background
Grafico FFT del power spectrum de la señal
Mediciones (por el momento): Vmax, Vmin, VPP, Periodo, Frecuencia fundamental (a futuro, THD, etc)
Cursor con mouse que informa VPP, Periodo y Frecuencia del area seleccionada en pantalla.

En fin, estoy muy contento. Lo que andaba haciendo ahora es ver la forma de pasar esto a esquematicos, no conozco mucho de programas para esto y hacerlo en livewire seria un infierno, sumado a que no existen todos los integrados ni en ese ni en proteus. Empece a hacer una especie de esquematico con photoshop, pegando las imagenes de los integrados en los datasheets, para unir con lineas pero ahi esta, pendiente.

El PCB esta adjuntado mas arriba aunque ya ha tenido cambios de componentes y algunas cirugias.

Adjunto aqui tambien el soft, tanto de la PC como el del PIC por si lo quieren chusmear o sacar ideas, o usarlo de base para algun hardware propio.

Esta hecho en vb.NET con framework 4. La comunicacion usb es via HID. Es un soft un tanto complejo y no esta muy comentado para entenderlo facil, pero cualquier duda, aqui estoy para lo que necesiten.
 

Adjuntos

  • Osciloscopio2013.zip
    1.3 MB · Visitas: 21
Atrás
Arriba