java - FFT 返回变为 NaN 的大值

Question

我正在使用 FFT 类来获取基频。我正在传递一些双值的数组。数组就像队列。当添加新值时数组将被更新。但我的问题是输出数组会不时变得很大。其变为E的幂值并最终返回NaN。我使用下面的 FFT 类，我很困惑问题出在哪里。如果任何人都可以通过找出原因来提供帮助，那将是一个很大的帮助。

这是我的 FFT 类(class)

  public class FFT {

  int n, m;

  // Lookup tables.  Only need to recompute when size of FFT changes.
  double[] cos;
  double[] sin;

  double[] window;

  public FFT(int n) {
    this.n = n;
    this.m = (int)(Math.log(n) / Math.log(2));

    // Make sure n is a power of 2
    if(n != (1<<m))
      throw new RuntimeException("FFT length must be power of 2");

    // precompute tables
    cos = new double[n/2];
    sin = new double[n/2];

//     for(int i=0; i<n/4; i++) {
//       cos[i] = Math.cos(-2*Math.PI*i/n);
//       sin[n/4-i] = cos[i];
//       cos[n/2-i] = -cos[i];
//       sin[n/4+i] = cos[i];
//       cos[n/2+i] = -cos[i];
//       sin[n*3/4-i] = -cos[i];
//       cos[n-i]   = cos[i];
//       sin[n*3/4+i] = -cos[i];    
//     }

    for(int i=0; i<n/2; i++) {
      cos[i] = Math.cos(-2*Math.PI*i/n);
      sin[i] = Math.sin(-2*Math.PI*i/n);
    }

    makeWindow();
  }

  protected void makeWindow() {
    // Make a blackman window:
    // w(n)=0.42-0.5cos{(2*PI*n)/(N-1)}+0.08cos{(4*PI*n)/(N-1)};
    window = new double[n];
    for(int i = 0; i < window.length; i++)
      window[i] = 0.42 - 0.5 * Math.cos(2*Math.PI*i/(n-1)) 
    + 0.08 * Math.cos(4*Math.PI*i/(n-1));
  }

  public double[] getWindow() {
    return window;
  }


  /***************************************************************
  * fft.c
  * Douglas L. Jones 
  * University of Illinois at Urbana-Champaign 
  * January 19, 1992 
  * http://cnx.rice.edu/content/m12016/latest/
  * 
  *   fft: in-place radix-2 DIT DFT of a complex input 
  * 
  *   input: 
  * n: length of FFT: must be a power of two 
  * m: n = 2**m 
  *   input/output 
  * x: double array of length n with real part of data 
  * y: double array of length n with imag part of data 
  * 
  *   Permission to copy and use this program is granted 
  *   as long as this header is included. 
  ****************************************************************/
  public void fft(double[] x, double[] y)
  {
    int i,j,k,n1,n2,a;
    double c,s,e,t1,t2;


    // Bit-reverse
    j = 0;
    n2 = n/2;
    for (i=1; i < n - 1; i++) {
      n1 = n2;
      while ( j >= n1 ) {
    j = j - n1;
    n1 = n1/2;
      }
      j = j + n1;

      if (i < j) {
    t1 = x[i];
    x[i] = x[j];
    x[j] = t1;
    t1 = y[i];
    y[i] = y[j];
    y[j] = t1;
      }
    }

    // FFT
    n1 = 0;
    n2 = 1;

    for (i=0; i < m; i++) {
      n1 = n2;
      n2 = n2 + n2;
      a = 0;

      for (j=0; j < n1; j++) {
    c = cos[a];
    s = sin[a];
    a +=  1 << (m-i-1);

    for (k=j; k < n; k=k+n2) {
      t1 = c*x[k+n1] - s*y[k+n1];
      t2 = s*x[k+n1] + c*y[k+n1];
      x[k+n1] = x[k] - t1;
      y[k+n1] = y[k] - t2;
      x[k] = x[k] + t1;
      y[k] = y[k] + t2;
    }
      }
    }
  }                          




  // Test the FFT to make sure it's working
  public static void main(String[] args) {
    int N = 8;

    FFT fft = new FFT(N);

    double[] window = fft.getWindow();
    double[] re = new double[N];
    double[] im = new double[N];

    // Impulse
    re[0] = 1; im[0] = 0;
    for(int i=1; i<N; i++)
      re[i] = im[i] = 0;
    beforeAfter(fft, re, im);

    // Nyquist
    for(int i=0; i<N; i++) {
      re[i] = Math.pow(-1, i);
      im[i] = 0;
    }
    beforeAfter(fft, re, im);

    // Single sin
    for(int i=0; i<N; i++) {
      re[i] = Math.cos(2*Math.PI*i / N);
      im[i] = 0;
    }
    beforeAfter(fft, re, im);

    // Ramp
    for(int i=0; i<N; i++) {
      re[i] = i;
      im[i] = 0;
    }
    beforeAfter(fft, re, im);

    long time = System.currentTimeMillis();
    double iter = 30000;
    for(int i=0; i<iter; i++)
      fft.fft(re,im);
    time = System.currentTimeMillis() - time;
    System.out.println("Averaged " + (time/iter) + "ms per iteration");
  }

  protected static void beforeAfter(FFT fft, double[] re, double[] im) {
    System.out.println("Before: ");
    printReIm(re, im);
    fft.fft(re, im);
    System.out.println("After: ");
    printReIm(re, im);
  }

  protected static void printReIm(double[] re, double[] im) {
    System.out.print("Re: [");
    for(int i=0; i<re.length; i++)
      System.out.print(((int)(re[i]*1000)/1000.0) + " ");

    System.out.print("]\nIm: [");
    for(int i=0; i<im.length; i++)
      System.out.print(((int)(im[i]*1000)/1000.0) + " ");

    System.out.println("]");
  }
}

下面是我在 android 中使用 FFT 实例的主要 Activity 类

    public class MainActivity extends Activity implements SensorEventListener{

    static final float ALPHA = 0.15f;

    private int count=0;
    private static GraphicalView view;
    private LineGraph line = new LineGraph();
    private static Thread thread;
    private SensorManager mSensorManager;
    private Sensor mAccelerometer;
    TextView title,tv,tv1,tv2,tv3,tv4,tv5,tv6;
    RelativeLayout layout;
    private double a;
    private double m = 0;
    private float p,q,r;
    public long[] myList;
    public double[] myList2;
    public double[] gettedList;
    static String k1,k2,k3,k4;
    int iniX=0;  
    public  FFT fft;
    public  myArray myArrayQueue;

    @Override
    protected void onCreate(Bundle savedInstanceState) {
        super.onCreate(savedInstanceState);
        setContentView(R.layout.activity_main);
         fft=new FFT(128);
         myList=new long[128];
         myList2=new double[128];
         gettedList=new double[128];
        myArrayQueue=new myArray(128);

        //get the sensor service
        mSensorManager = (SensorManager) getSystemService(Context.SENSOR_SERVICE);
        //get the accelerometer sensor
        mAccelerometer = mSensorManager.getDefaultSensor(Sensor.TYPE_ACCELEROMETER);
        //get layout
        layout = (RelativeLayout)findViewById(R.id.relative);
        LinearLayout layout = (LinearLayout) findViewById(R.id.layoutC);
        view= line.getView(this);
        layout.addView(view);
        //get textviews
        title=(TextView)findViewById(R.id.name);
        //tv=(TextView)findViewById(R.id.xval);
        //tv1=(TextView)findViewById(R.id.yval);
        //tv2=(TextView)findViewById(R.id.zval);


        tv3=(TextView)findViewById(R.id.TextView04);
        tv4=(TextView)findViewById(R.id.TextView01);
        tv5=(TextView)findViewById(R.id.TextView02);
        tv6=(TextView)findViewById(R.id.TextView03);



        for (int i = 0; i < myList2.length; i++){
            myList2[i] =0;

        }



    }

    public final void onAccuracyChanged(Sensor sensor, int accuracy)
    {
        // Do something here if sensor accuracy changes.
    }
    @Override
    public final void onSensorChanged(SensorEvent event)
    {
        count=+1;
        // Many sensors return 3 values, one for each axis.


        float x = event.values[0];
        float y = event.values[1];
        float z = event.values[2];



        //float[] first={x,y,z};
        //  float[] larst={p,q,r};

        //larst= lowPass(first,larst);
        //double FY= b.Filter(y);
        //double FZ= b.Filter(z);

        //get merged value
        //  m = (float) Math.sqrt(larst[0]*larst[0]+larst[1]*larst[1]+larst[2]*larst[2]);
        m=(double)Math.sqrt(x*x+y*y+z*z);


        //display values using TextView
        //title.setText(R.string.app_name);
        //tv.setText("X axis" +"\t\t"+x);
        //tv1.setText("Y axis" + "\t\t" +y);
        //tv2.setText("Z axis" +"\t\t" +z);



        //myList[iniX]=m*m;
        //myList[iniX+1]=myList[iniX];

            iniX=+1;

            //myList[3]=myList[2];
            //myList[2]=myList[1];
            //myList[1]=myList[0];

        myArrayQueue.insert(m*m);
    gettedList=myArrayQueue.getMyList();
        /*  for(int a = myList.length-1;a>0;a--)
            {
            myList[a]=myList[a-1];


            }
            myList[0]=m*m;
        */  


    fft.fft(gettedList, myList2);
            k1=Double.toString(myList2[0]);
            k2=Double.toString(myList2[1]);
            k3=Double.toString(myList2[2]);
            k4=Double.toString(myList2[3]);

            tv3.setText("[0]= "+k1);
            tv4.setText("[1]= "+k2);
            tv5.setText("[2]= "+k3);
            tv6.setText("[3]= "+k4);
            line.addNewPoint(iniX,(float) m);
            view.repaint();


    }

    @Override
    protected void onResume()
    {
        super.onResume();
        mSensorManager.registerListener(this, mAccelerometer, SensorManager.SENSOR_DELAY_NORMAL);
    }
    @Override
    protected void onPause()
    {
        super.onPause();
        mSensorManager.unregisterListener(this);
    }

    public void LineGraphHandler(View view){



    }

    //Low pass filter
    protected float[] lowPass( float[] input, float[] output ) {
        if ( output == null ) return input;

        for ( int i=0; i<input.length; i++ ) {
            output[i] = output[i] + ALPHA * (input[i] - output[i]);
        }
        return output;
    }
    /*@Override
    public void onStart(){
        super.onStart();    
        view= line.getView(this);
        setContentView(view);
    }*/

}

Answer 1

如果输入包含 NaN，FFT 输出只会产生 NaN。因此，在调用 FFT 来调试代码之前，请显式检查 FFT 的输入数组是否有任何超出范围的值。然后从那里倒推找出它们来自哪里。