Estoy intentando implementar un filtro de parada de banda (muesca) en un microcontrolador dsPIC33EP64GS506 para filtrar un componente 100 Hz de una señal de entrada. El problema aquí es que mi microcontrolador no tiene una unidad de punto flotante, por lo que tengo que usar una aritmética de punto fijo, pero para ser más precisos, uso la aritmética de enteros. El tiempo de muestra del sistema es \ $ T_s = 50 ~ \ mu \ text {s} \ $. Aquí explico mi problema en detalle, y las preguntas están al final del post. Encontrará adjunto un código MATLAB para ejecutar simulaciones si desea: enlace de descarga en mi Dropbox . Tenga en cuenta que no tiene que ser un usuario registrado de Dropbox para poder descargar el archivo.
La resolución ADC es 12b , que "aumente" a 15b , no para aumentar la resolución de medición (que no se puede hacer), sino para aumentar la resolución del filtro:
v_in = ADCBUF0<<3;
La función de transferencia de filtro de muesca en el dominio s se da como:
$$ G_s (s) = \ frac {s ^ 2 + 2 \ zeta_1 \ omega_n s + \ omega_n ^ 2} {s ^ 2 + 2 \ zeta_2 \ omega_ns + \ omega_n ^ 2}, $$
donde los parámetros de filtro \ $ \ zeta_1 \ $, \ $ \ zeta_2 \ $, y \ $ \ omega_n \ $ determinan completamente el filtro. Si queremos filtrar un componente 100 Hz , entonces configuramos \ $ \ omega_n = 2 \ pi \ cdot100 ~ \ text {s} ^ {- 1} \ $. En cuanto a otros parámetros, sin explicaciones detalladas, uso \ $ \ zeta_1 = 0.001 \ $ y \ $ \ zeta_2 = 1 \ $. Para implementar este filtro en un sistema digital, necesitamos discretizar la función de transferencia en s-domain, y para eso, uso un Método de discretización de Tustin .
$$ G_z (z) = \ frac {b_0 + b_1 z ^ {- 1} + b_2 z ^ {- 2}} {1 + a_1 z ^ {- 1} + a_2 z ^ {- 2}} . $$
La ecuación recursiva correspondiente, tal como se implementa en un sistema digital, es:
$$ y (k) = b_0 u (k) + b_1 u (k-1) + b_2 u (k-2) - a_1 y (k-1) - a_2 y (k-2). $$
Aquí está el código MATLAB para obtener la función de transferencia en el dominio z:
Ts = 50e-6;
zeta1 = 0.001;
zeta2 = 1;
wn = 100*(2*pi);
Gs = tf([1 2*zeta1*wn wn^2], [1 2*zeta2*wn wn^2]);
Gz = c2d(Gs, Ts, 'tustin');
bodeplot(Gz);
Consulte a continuación las características de frecuencia del filtro en el dominio z generado por MATLAB. Como se puede ver, este tipo de filtro es muy sensible en términos de una respuesta de frecuencia. Incluso los cambios más pequeños en los parámetros podrían causar un comportamiento completamente diferente, por ejemplo, una ganancia inesperada, un cambio de fase, etc.
Ahora,comonotengounaunidaddepuntoflotantedisponible,utilizounmétodoconocidollamadoescaladobinario.Cualquiernúmerodecimal\$d\in\mathbb{R}\$puederepresentarsecomo\$\frac{\lfloord\cdot2^r\rceil}{2^r}\$,donde\$r\in\mathbb{N}\$,y\$\lfloor\cdot\rceil\$esunafuncióndeenteroredondomáscercano.Elescaladobinarioesunmétodoutilizadocuandoqueremosevitarelusodeunadivisión,queesunaoperaciónmuycostosaentérminosdeciclosdeCPUrequeridos(típicamente18ciclos),yaqueeldesplazamientobinariopuedeproducirlosmismosresultadosenmuchosmenosciclos(típicamente1ciclo).Paraesteejemplousé\$r=15\$,queeslaprecisiónmásaltaquepodríausarconsiderandolosbitsdisponiblesparahacerunamultiplicación(unresultadodebeestardentrodelos32bits).Laecuaciónrecursivacorrespondientequeutilizaaritméticadeenterosseimplementadelasiguientemanera:
yk0=(B0*uk0+B1*uk1+B2*uk2-A1*yk1-A2*yk2+(1<<14))>>15;uk2=uk1;uk1=uk0;yk2=yk1;yk1=yk0;Eltérmino(1<<14)seusaparaasegurarqueelresultadoseredondeaalenteromáscercanodespuésdelaoperacióndecambiodebit.Tengaencuentaque>>15esprácticamenteunadivisiónenterapor32768,perosabemosquenoescompletamenteequivalente.Eldesplazamientodebitssiempreseredondeaamenosinfinito,mientrasqueladivisiónenterasiempreseredondeaacero.Losvaloresdelosparámetrosdelfiltrosonlossiguientes:B0=31771,B1=-63509mB2=31769,A1=-63509,A2=30772.EstotambiénsepuedeencontrarenuncódigoMATLABproporcionado.
Mesorprendiómuchocuandomedicuentadequelavariantedeenterosdelfiltronofuncionaenabsoluto.Veaacontinuaciónlasrespuestasdetresfiltrosdiferentesutilizadosparaeliminaruncomponente100Hzdelaseñaldeentrada.Dearribaaabajo:(1)filtrodemuescaenlaimplementacióndepuntoflotante,(2)filtrodemuescaenlaimplementaciónaritméticadeenterosusando\$r=15\$,(3)unfiltrodepromediomóvilsimple.
Aquí es cómo he implementado un filtro de promedio móvil en C:
uk200 = window[ind];
window[ind] = uk0;
win_sum = win_sum - uk200 + uk0;
yk0 = (((win_sum+(1<<3))>>4)*5243+(1<<15))>>16;
ind++;
if (ind==200) ind=0;
Preguntas para la comunidad.
Como no tengo mucha experiencia en el filtrado digital, ¿puede alguien confirmarme estos resultados esperados? ¿Es realmente tan problemático implementar un filtro de muesca en un sistema digital usando solo aritmética de enteros ? ¿Qué filtro digital suele utilizarse para eliminar una determinada frecuencia ? Como veo en estos resultados, el filtro de promedio móvil supera a ambas implementaciones de filtros de muesca. ¿Hay algún inconveniente para el filtro de promedio móvil que debería tener en cuenta, a excepción del aumento evidente de la demanda de memoria ?