yixiaofan 发表于 2007-4-12 22:43


明白了,谢谢

zhaoyixu 发表于 2007-6-18 11:09

现在我是一点也看不懂,不过我相信以后就会懂了。

xiaodinger 发表于 2008-2-27 10:33

hr(1)=(wl-wh)/pi;
hr(2:M+1)=(sin(wl*k)-sin(wh*k))./(pi*k).*w;
hi(1)=0;
hi(2:M+1)=(cos(wl*k)-cos(wh*k))./(pi*k).*w;

程序是不是有误,对照丁康的论文,缺少了k=1时的值,却多了k=M+1时的值?

xiaodinger 发表于 2008-2-27 10:56

晕了,看明白了,不好意思

晕了,看明白了,不好意思

JLBhaidao 发表于 2009-6-17 03:58

zoomfft函数不能调用

怎么会出现“??? Undefined command/function 'zoomfft'.“,高手教我!

ChaChing 发表于 2009-6-17 08:17

回复 50楼 JLBhaidao 的帖子

使用which zoomfft检查下
不是无此函数就是路径没设好

diandidiandi 发表于 2010-2-10 20:10

弱弱的问一下
for k=1:N
      kk=(k-1)*D+M+(i-1)*N;
      xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
      xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(x(kk+2:kk+M+1)-x(kk:-1:kk-M+1)));
    end
这段程序是在进行滤波和抽取吗?

无限剑制 发表于 2010-2-11 14:17

学习了:lol

zmp007 发表于 2010-4-9 22:30

回复 楼主 yangzj 的帖子

请教给位仁兄,我运行这个程序时,会出现以下错误信号,不知为何?请指教,多谢!
??? function =ZoomFFT(x,fs,N,fe,D,L,M)
    |
Error: Function definitions are not permitted at the prompt or in scripts.

zxk5511 发表于 2010-8-4 21:05

学习了:@L 还是不太懂

secondye 发表于 2010-8-24 00:07

谢谢楼主,拷下来学习了

土木年华 发表于 2010-8-26 11:30

回复 shanghai 的帖子


   

200810100659 发表于 2010-11-30 16:44

{:{39}:}

ylhlly 发表于 2011-4-16 13:18

纠结这个细化谱好几天了,也找了网上和书上其它的matlab源程序,但改写C#始终不对。。。今天用楼主的这个就对了~~非常感谢。。

ylhlly 发表于 2011-4-16 13:23

    private ArrayList m_TimeSeries = new ArrayList();         //存储时域信号的数据序列
    private ArrayList m_ZoomFFTSeries = new ArrayList();      //存储细化谱的数据序列
    //tempTr、tempTi暂存kbfft()的入口数据;tempFr、tempFi存储kbfft()返回的实部和虚部
    private ArrayList tempTr = new ArrayList(), tempTi = new ArrayList(), tempFr = new ArrayList(), tempFi = new ArrayList();
    private double m_Xstep, m_ZoomXmax, m_ZoomXmin;

    public ArrayList TimeSeries
    {
      get{ return m_TimeSeries; }
      set{ m_TimeSeries = value; }
    }
    public ArrayList ZoomFFTSeries
    {
      get { return m_ZoomFFTSeries; }
      set { m_ZoomFFTSeries = value; }
    }
    public double Xstep
    {
      get { return m_Xstep; }
      set { m_Xstep = value; }
    }
    public double ZoomXmax
    {
      get { return m_ZoomXmax; }
      set { m_ZoomXmax = value; }
    }
    public double ZoomXmin
    {
      get { return m_ZoomXmin; }
      set { m_ZoomXmin = value; }
    }


    #region 细化谱序列
    /// <summary>
    /// 通过时域信号的数据序列,计算细化谱数据序列。fin为分析中心频率,iN为做谱点数,scale细化倍数,iL平均段数,iM为滤波器半阶数
    /// </summary>
    public void GetZoomFFTSeries(double fin, int scale, int iN, int iL, int iM)//
    {
      m_ZoomFFTSeries.Clear();
      double[] w = new double;
      for (int i = 0; i < iM; i++)
      {
            w = 0.5 + 0.5 * Math.Cos(Math.PI * (i+1) / iM);//Hanning窗
      }
      double fl, fh;//最小截止频率、最大截止频率
      //fl = (fin - m_SampleHz / (4.0 * scale)) > (-m_SampleHz / 2.2) ? (fin - m_SampleHz / (4 * scale)) : -m_SampleHz / 2.2;
      //fh = (fin + m_SampleHz / (4.0 * scale)) < (m_SampleHz / 2.2) ? (fin + m_SampleHz / (4 * scale)) : m_SampleHz / 2.2;
      fl = fin - m_SampleHz / (4.0 * scale); fh = fin + m_SampleHz / (4.0 * scale);

      double yf = scale * fl;//移频量
      double df = m_SampleHz / scale / iN;//

      double[] xz = new double;
      double wl = 2 * Math.PI * fl / m_SampleHz, wh = 2 * Math.PI * fh / m_SampleHz;
      double[] hr = new double, hi = new double;
      hr = (wl - wh) / Math.PI; hi = 0;
      for (int i = 1; i <= iM; i++)
      {
            hr = (Math.Sin(wl * i) - Math.Sin(wh * i)) / (Math.PI * i) * w;
            hi = (Math.Cos(wl * i) - Math.Cos(wh * i)) / (Math.PI * i) * w;
      }
      double[] w1 = new double;
      for (int i = 0; i < iN; i++)
      {
            w1 = 0.5 - 0.5 * Math.Cos(2 * Math.PI * i / iN);
      }

      double[] xrz = new double, xiz = new double,xzt=new double;
      double[] x = new double;
      for (int i = 0; i < m_TimeSeries.Count; i++)
      {
            x = Convert.ToDouble(m_TimeSeries);
      }
      double factor = -2 * Math.PI * yf / m_SampleHz;
      for (int i = 1; i <= iL; i++)
      {
            for (int k = 1; k <= iN; k++)
            {
                int kk = (k - 1) * scale + iM + (i - 1) * iN;
                double sumR = 0, sumI = 0;
                for (int j = 1; j <= iM; j++)
                {
                  sumR += hr * (x + x);
                  sumI += hi * (x - x);
                }
                xrz=(x * hr + sumR);
                xiz=(x * hi + sumI);
            }
            double sumXZTr = 0, sumXZTi = 0;
            for (int k = 0; k < 2*iN; k+=2)
            {
                xzt = xrz * Math.Cos(k * factor / 2) - xiz * Math.Sin(k * factor / 2);//偶数位为实部,奇数位为虚部
                xzt = xiz * Math.Cos(k * factor / 2) + xrz * Math.Sin(k * factor / 2);
                xzt = xzt * w1;
                xzt = xzt * w1;
                sumXZTr += xzt; sumXZTi += xzt;
            }
            tempTr.Clear(); tempTi.Clear();
            for (int k = 0; k < 2 * iN; k += 2)
            {
                tempTr.Add(xzt - sumXZTr / iN);
                tempTi.Add(xzt - sumXZTi / iN);
            }
            kbfft(tempTr.Count, 0, 0);
            for (int k = 0; k < iN / 2; k++)
            {
                double temR = Convert.ToDouble(tempFr), temI = Convert.ToDouble(tempFi);
                xz += Math.Pow(Math.Sqrt(temI * temI + temR * temR) / iN * 2, 2);
            }
      }
      for (int i = 0; i < iN / 2; i++)
      {
            m_ZoomFFTSeries.Add(Math.Sqrt(xz / iL));
      }
      m_ZoomXmin = fl; m_ZoomXmax = fh;
      m_Xstep = (fh - fl) / iN * 2;
    }
   
    /// <summary>
    /// GetZoomFFTSeries(double fin, int scale, int iN, int iL, int iM)的重构方法。
    /// </summary>
    public void GetZoomFFTSeries()
    {
      this.GetZoomFFTSeries(500, 5, 1024, 1, m_TimeSeries.Count/16);
    }
    #endregion

    #region 快速傅里叶变换方法
    /***************************************************************************
   * 入口参数:
   *l: l = 0, 傅立叶变换; l = 1, 逆傅立叶变换
   * il: il = 0,不计算傅立叶变换或逆变换模和幅角;il = 1,计算模和幅角
   * n: 输入的点数,为偶数,一般为32,64,128,...,1024等
   * k: 满足n=2^k(k>0),实质上k是n个采样数据可以分解为偶次幂和奇次幂的次数
   * pr[]: l=0时,存放N点采样数据的实部
   * l=1时, 存放傅立叶变换的N个实部
   * pi[]: l=0时,存放N点采样数据的虚部
   * l=1时, 存放傅立叶变换的N个虚部
   *
   * 出口参数:
   * fr[]: l=0, 返回傅立叶变换的实部
   * l=1, 返回逆傅立叶变换的实部
   * fi[]: l=0, 返回傅立叶变换的虚部
   * l=1, 返回逆傅立叶变换的虚部
   * pr[]: il = 1,l = 0 时,返回傅立叶变换的模
   * il = 1,l = 1 时,返回逆傅立叶变换的模
   * pi[]: il = 1,l = 0 时,返回傅立叶变换的辐角
   * il = 1,l = 1 时,返回逆傅立叶变换的辐角
   ***************************************************************************/
    /// <summary>
    /// 快速傅里叶变换方法。
    /// </summary>
    void kbfft(int n,int l, int il)
    {
      int k = (int)Math.Log(n, 2);
      double[] pr = new double, pi = new double, fr = new double, fi = new double;
      for (int i = 0; i < n; i++)
      {
            if (i < tempTr.Count)
                pr = Convert.ToDouble(tempTr);
            else
                pr = 0;
            if (i < tempTi.Count)
                pi = Convert.ToDouble(tempTi);
            else
                pi = 0;
      }

      intm, is1, j, nv;
      double p, q, s, vr, vi, poddr, poddi;

      //排序
      for (int it = 0; it <= n - 1; it++)
      {
            m = it; is1 = 0;
            for (int i = 0; i <= k - 1; i++)
            {
                j = m / 2; is1 = 2 * is1 + (m - 2 * j); m = j;
                fr = pr; fi = pi;
            }
      }
      
      //蝶形运算
      pr = 1.0; pi = 0.0;
      p = Math.PI * 2 / (1.0 * n);
      pr = Math.Cos(p); pi = -Math.Sin(p);
      if (l != 0) pi = -pi;
      for (int i = 2; i <= n - 1; i++)
      {
            p = pr * pr; q = pi * pi;
            s = (pr + pi) * (pr + pi);
            pr = p - q; pi = s - p - q;
      }
      for (int it = 0; it <= n - 2; it = it + 2)
      {
            vr = fr; vi = fi;
            fr = vr + fr; fi = vi + fi;
            fr = vr - fr; fi = vi - fi;
      }
      m = n / 2; nv = 2;
      for (int l0 = k - 2; l0 >= 0; l0--)
      {
            m = m / 2; nv = 2 * nv;
            for (int it = 0; it <= (m - 1) * nv; it = it + nv)
                for ( j = 0; j <= (nv / 2) - 1; j++)
                {
                  p = pr * fr;
                  q = pi * fi;
                  s = pr + pi;
                  s = s * (fr + fi);
                  poddr = p - q; poddi = s - p - q;
                  fr = fr - poddr;
                  fi = fi - poddi;
                  fr = fr + poddr;
                  fi = fi + poddi;
                }
      }
      if (l != 0)
            for (int i = 0; i <= n - 1; i++)
            {
                fr = fr / (1.0 * n);
                fi = fi / (1.0 * n);
            }
      if (il != 0)
            for (int i = 0; i <= n - 1; i++)
            {
                pr = Math.Sqrt(fr * fr + fi * fi);
                pr = (pr / (n / 2)); //各次谐波幅值,其中pr为基波幅值
                if (Math.Abs(fr) < 0.000001 * Math.Abs(fi))
                {
                  if ((fi * fr) > 0) pi = 90.0;
                  else pi = -90.0;
                }
                else
                  pi = Math.Atan(fi / fr) * 360.0 / (2*Math.PI);
            }


      tempFi.Clear(); tempFr.Clear();
      for (int i = 0; i < n; i++)
      {
            tempFr.Add(fr);
            tempFi.Add(fi);
      }
    }
    #endregion这是我写的一个类(其中一部分)。。。希望对大家有所帮助
页: 1 2 3 [4] 5
查看完整版本: [原创]zoomfft