java - studio - programacion android pdf 2018
Cómo implementar el filtro de paso bajo usando java. (6)
Estoy tratando de implementar un filtro de paso bajo en Java. Mi requerimiento es muy simple, tengo que eliminar las señales más allá de una frecuencia particular (dimensión única). Parece que el filtro Butterworth se adapta a mi necesidad.
Ahora lo importante es que el tiempo de CPU debe ser lo más bajo posible. Habría cerca de un millón de muestras que el filtro tendría que procesar y a nuestros usuarios no les gusta esperar demasiado. ¿Hay alguna implementación prefabricada de filtros Butterworth que tenga algoritmos óptimos para el filtrado?
Adopté esto de http://www.dspguide.com/ Soy bastante nuevo en java, así que no es bonito, pero funciona
/*
* To change this license header, choose License Headers in Project Properties.
* To change this template file, choose Tools | Templates
* and open the template in the editor.
*/
package SoundCruncher;
import java.util.ArrayList;
/**
*
* @author 2sloth
* filter routine from "The scientist and engineer''s guide to DSP" Chapter 20
* filterOrder can be any even number between 2 & 20
* cutoffFreq must be smaller than half the samplerate
* filterType: 0=lowPass 1=highPass
* ripplePercent is amount of ripple in Chebyshev filter (0-29) (0=butterworth)
*/
public class Filtering {
double[] filterSignal(ArrayList<Float> signal, double sampleRate ,double cutoffFreq, double filterOrder, int filterType, double ripplePercent) {
double[][] recursionCoefficients = new double[22][2];
// Generate double array for ease of coding
double[] unfilteredSignal = new double[signal.size()];
for (int i=0; i<signal.size(); i++) {
unfilteredSignal[i] = signal.get(i);
}
double cutoffFraction = cutoffFreq/sampleRate; // convert cut-off frequency to fraction of sample rate
System.out.println("Filtering: cutoffFraction: " + cutoffFraction);
//ButterworthFilter(0.4,6,ButterworthFilter.Type highPass);
double[] coeffA = new double[22]; //a coeffs
double[] coeffB = new double[22]; //b coeffs
double[] tA = new double[22];
double[] tB = new double[22];
coeffA[2] = 1;
coeffB[2] = 1;
// calling subroutine
for (int i=1; i<filterOrder/2; i++) {
double[] filterParameters = MakeFilterParameters(cutoffFraction, filterType, ripplePercent, filterOrder, i);
for (int j=0; j<coeffA.length; j++){
tA[j] = coeffA[j];
tB[j] = coeffB[j];
}
for (int j=2; j<coeffA.length; j++){
coeffA[j] = filterParameters[0]*tA[j]+filterParameters[1]*tA[j-1]+filterParameters[2]*tA[j-2];
coeffB[j] = tB[j]-filterParameters[3]*tB[j-1]-filterParameters[4]*tB[j-2];
}
}
coeffB[2] = 0;
for (int i=0; i<20; i++){
coeffA[i] = coeffA[i+2];
coeffB[i] = -coeffB[i+2];
}
// adjusting coeffA and coeffB for high/low pass filter
double sA = 0;
double sB = 0;
for (int i=0; i<20; i++){
if (filterType==0) sA = sA+coeffA[i];
if (filterType==0) sB = sB+coeffB[i];
if (filterType==1) sA = sA+coeffA[i]*Math.pow(-1,i);
if (filterType==1) sB = sB+coeffA[i]*Math.pow(-1,i);
}
// applying gain
double gain = sA/(1-sB);
for (int i=0; i<20; i++){
coeffA[i] = coeffA[i]/gain;
}
for (int i=0; i<22; i++){
recursionCoefficients[i][0] = coeffA[i];
recursionCoefficients[i][1] = coeffB[i];
}
double[] filteredSignal = new double[signal.size()];
double filterSampleA = 0;
double filterSampleB = 0;
// loop for applying recursive filter
for (int i= (int) Math.round(filterOrder); i<signal.size(); i++){
for(int j=0; j<filterOrder+1; j++) {
filterSampleA = filterSampleA+coeffA[j]*unfilteredSignal[i-j];
}
for(int j=1; j<filterOrder+1; j++) {
filterSampleB = filterSampleB+coeffB[j]*filteredSignal[i-j];
}
filteredSignal[i] = filterSampleA+filterSampleB;
filterSampleA = 0;
filterSampleB = 0;
}
return filteredSignal;
}
/* pi=3.14...
cutoffFreq=fraction of samplerate, default 0.4 FC
filterType: 0=LowPass 1=HighPass LH
rippleP=ripple procent 0-29 PR
iterateOver=1 to poles/2 P%
*/
// subroutine called from "filterSignal" method
double[] MakeFilterParameters(double cutoffFraction, int filterType, double rippleP, double numberOfPoles, int iteration) {
double rp = -Math.cos(Math.PI/(numberOfPoles*2)+(iteration-1)*(Math.PI/numberOfPoles));
double ip = Math.sin(Math.PI/(numberOfPoles*2)+(iteration-1)*Math.PI/numberOfPoles);
System.out.println("MakeFilterParameters: ripplP:");
System.out.println("cutoffFraction filterType rippleP numberOfPoles iteration");
System.out.println(cutoffFraction + " " + filterType + " " + rippleP + " " + numberOfPoles + " " + iteration);
if (rippleP != 0){
double es = Math.sqrt(Math.pow(100/(100-rippleP),2)-1);
// double vx1 = 1/numberOfPoles;
// double vx2 = 1/Math.pow(es,2)+1;
// double vx3 = (1/es)+Math.sqrt(vx2);
// System.out.println("VX''s: ");
// System.out.println(vx1 + " " + vx2 + " " + vx3);
// double vx = vx1*Math.log(vx3);
double vx = (1/numberOfPoles)*Math.log((1/es)+Math.sqrt((1/Math.pow(es,2))+1));
double kx = (1/numberOfPoles)*Math.log((1/es)+Math.sqrt((1/Math.pow(es,2))-1));
kx = (Math.exp(kx)+Math.exp(-kx))/2;
rp = rp*((Math.exp(vx)-Math.exp(-vx))/2)/kx;
ip = ip*((Math.exp(vx)+Math.exp(-vx))/2)/kx;
System.out.println("MakeFilterParameters (rippleP!=0):");
System.out.println("es vx kx rp ip");
System.out.println(es + " " + vx*100 + " " + kx + " " + rp + " " + ip);
}
double t = 2*Math.tan(0.5);
double w = 2*Math.PI*cutoffFraction;
double m = Math.pow(rp, 2)+Math.pow(ip,2);
double d = 4-4*rp*t+m*Math.pow(t,2);
double x0 = Math.pow(t,2)/d;
double x1 = 2*Math.pow(t,2)/d;
double x2 = Math.pow(t,2)/d;
double y1 = (8-2*m*Math.pow(t,2))/d;
double y2 = (-4-4*rp*t-m*Math.pow(t,2))/d;
double k = 0;
if (filterType==1) {
k = -Math.cos(w/2+0.5)/Math.cos(w/2-0.5);
}
if (filterType==0) {
k = -Math.sin(0.5-w/2)/Math.sin(w/2+0.5);
}
d = 1+y1*k-y2*Math.pow(k,2);
double[] filterParameters = new double[5];
filterParameters[0] = (x0-x1*k+x2*Math.pow(k,2))/d; //a0
filterParameters[1] = (-2*x0*k+x1+x1*Math.pow(k,2)-2*x2*k)/d; //a1
filterParameters[2] = (x0*Math.pow(k,2)-x1*k+x2)/d; //a2
filterParameters[3] = (2*k+y1+y1*Math.pow(k,2)-2*y2*k)/d; //b1
filterParameters[4] = (-(Math.pow(k,2))-y1*k+y2)/d; //b2
if (filterType==1) {
filterParameters[1] = -filterParameters[1];
filterParameters[3] = -filterParameters[3];
}
// for (double number: filterParameters){
// System.out.println("MakeFilterParameters: " + number);
// }
return filterParameters;
}
}
Aquí hay un filtro de paso bajo que utiliza una transformación de Fourier en la biblioteca matemática de Apache.
public double[] fourierLowPassFilter(double[] data, double lowPass, double frequency){
//data: input data, must be spaced equally in time.
//lowPass: The cutoff frequency at which
//frequency: The frequency of the input data.
//The apache Fft (Fast Fourier Transform) accepts arrays that are powers of 2.
int minPowerOf2 = 1;
while(minPowerOf2 < data.length)
minPowerOf2 = 2 * minPowerOf2;
//pad with zeros
double[] padded = new double[minPowerOf2];
for(int i = 0; i < data.length; i++)
padded[i] = data[i];
FastFourierTransformer transformer = new FastFourierTransformer(DftNormalization.STANDARD);
Complex[] fourierTransform = transformer.transform(padded, TransformType.FORWARD);
//build the frequency domain array
double[] frequencyDomain = new double[fourierTransform.length];
for(int i = 0; i < frequencyDomain.length; i++)
frequencyDomain[i] = frequency * i / (double)fourierTransform.length;
//build the classifier array, 2s are kept and 0s do not pass the filter
double[] keepPoints = new double[frequencyDomain.length];
keepPoints[0] = 1;
for(int i = 1; i < frequencyDomain.length; i++){
if(frequencyDomain[i] < lowPass)
keepPoints[i] = 2;
else
keepPoints[i] = 0;
}
//filter the fft
for(int i = 0; i < fourierTransform.length; i++)
fourierTransform[i] = fourierTransform[i].multiply((double)keepPoints[i]);
//invert back to time domain
Complex[] reverseFourier = transformer.transform(fourierTransform, TransformType.INVERSE);
//get the real part of the reverse
double[] result = new double[data.length];
for(int i = 0; i< result.length; i++){
result[i] = reverseFourier[i].getReal();
}
return result;
}
Como dijo Mark Peters en su comentario: Un filtro que necesita filtrar mucho debe escribirse en C o C ++. Pero todavía puedes hacer uso de Java. Simplemente eche un vistazo a la interfaz nativa de Java (JNI) . Debido a que C / C ++ se compila al código de máquina nativo, se ejecutará mucho más rápido que ejecutar su código de bytes en la Máquina Virtual de Java (JVM), que de hecho es un procesador virtual que traduce el código de bytes a la máquina local su código nativo (dependiendo de en el conjunto de instrucciones de CPU como x86, x64, ARM , ....)
El diseño de filtros es un arte de compensaciones, y para hacerlo bien debe tener en cuenta algunos detalles.
¿Cuál es la frecuencia máxima que debe pasarse "sin mucha" atención, y cuál es el valor máximo de "sin mucha"?
¿Cuál es la frecuencia mínima que debe atenuarse "mucho" y cuál es el valor mínimo de "mucho"?
¿Cuánta ondulación (es decir, la variación en la atenuación) es aceptable dentro de las frecuencias que se supone que pasa el filtro?
Tiene una amplia gama de opciones, que le costarán una variedad de cantidades de cómputo. Un programa como matlab o scilab puede ayudarlo a comparar las compensaciones . Querrá familiarizarse con conceptos como la expresión de frecuencias como una fracción decimal de una frecuencia de muestreo, y el intercambio entre las mediciones de atenuación lineal y log (dB).
Por ejemplo, un filtro de paso bajo "perfecto" es rectangular en el dominio de la frecuencia. Expresado en el dominio del tiempo como una respuesta de impulso, sería una función sinc (sin x / x) con las colas alcanzando el infinito positivo y negativo. Obviamente, no puede calcular eso, por lo que la pregunta es si aproxima la función sinc a una duración finita que puede calcular, ¿cuánto de eso degradará su filtro?
Alternativamente, si desea un filtro de respuesta de impulso finito que sea muy barato de calcular, puede usar un "carro de caja" o un filtro rectangular donde todos los coeficientes son 1. (Esto puede hacerse incluso más barato si lo implementa como un filtro CIC explotando el desbordamiento binario para hacer acumuladores ''circulares'', ya que de todos modos tomará el derivado más adelante). Pero un filtro que es rectangular en el tiempo parece una función sinc en la frecuencia: tiene un rollo sin x / x en la banda de paso (a menudo elevado a cierta potencia, ya que normalmente tendría una versión de múltiples etapas), y algunos "rebotan" en la banda de parada. Aún en algunos casos es útil, ya sea por sí mismo o cuando es seguido por otro tipo de filtro.
Recientemente he diseñado una función simple de Butterworth ( http://baumdevblog.blogspot.com/2010/11/butterworth-lowpass-filter-coefficients.html ). Son fáciles de codificar en Java y deberían ser lo suficientemente rápidos si me preguntas (solo tendrías que cambiar el filtro (doble * muestras, int count) para filtrar (doble [] muestras, int count), supongo).
El problema con JNI es que cuesta independencia de la plataforma, puede confundir al compilador del punto de acceso y las llamadas al método JNI dentro de su código aún pueden ralentizar las cosas. Así que recomendaría probar Java y ver si es lo suficientemente rápido.
En algunos casos, puede ser beneficioso usar primero una transformada rápida de Fourier y aplicar el filtrado en el dominio de la frecuencia, pero dudo que esto sea más rápido que aproximadamente 6 multiplicados y algunas adiciones por muestra para un filtro de paso bajo simple.
Tengo una página que describe un filtro de paso bajo de CPU muy simple y muy bajo que también puede ser independiente de la tasa de cuadros. Lo uso para suavizar la entrada del usuario y también para graficar las tasas de cuadros a menudo.
http://phrogz.net/js/framerate-independent-low-pass-filter.html
En resumen, en su bucle de actualización:
// If you have a fixed frame rate
smoothedValue += (newValue - smoothedValue) / smoothing
// If you have a varying frame rate
smoothedValue += timeSinceLastUpdate * (newValue - smoothedValue) / smoothing
Un valor de smoothing de 1 hace que no se produzca un suavizado, mientras que valores más altos suavizan cada vez más el resultado.
La página tiene un par de funciones escritas en JavaScript, pero la fórmula es independiente del lenguaje.