在一条曲线上检出多个峰的问题,一直困扰着我,现在时间越来越紧了。还是没有能够相处一个解决方法。
请教算法,怎样检出曲线上的多个峰??敬请赐教,多谢了!

解决方案 »

  1.   

    也就是求极大值点
    在极大值点
    其导函数值为0
    且导函数在极大值附近是递减的如
    y = -x*x 其导函数为 y = -2*x当x = 0时,y=-2*0 = 0,导函数值为0,且y=-2x是递减的
    其极大值点就是(0, 0);
      

  2.   

    用三次样条插值法,根据已知的点求出平滑曲线,然后求导数值,导数值为零的地方即为波峰或波谷,到底是波峰还是波谷要根据该点的函数值与附近点的函数值比较而定,这一点很简单。附上三次样条插值C语言函数:/***********************************************************************************************
    // 
    //【函  数】double Espl1( double *x, double *y, 
    // int n, double *dy, 
    // double *ddy, double *t, 
    // int m, double *z, 
    // double *dz, double *ddz ) 
    //【参  数】x     : 双精度实型一维数组,长度为n,存放n个给定结点的值x(i);要求x(0)<x(1)<....<x(n-1)。
    //          y     : 双精度实型一维数组,长度为n,存放n个给定结点的函数值y(i),y(i)=f(x(i)),i=0,1,...,n-1,
    //          n     : 整形变量,给定结点的个数,        
    //          dy    : 双精度实型一维数组,长度为n,调用时,dy[0]存放给定的左端点处的一阶导数值,dy[n-1]存放
    // 给定的右端点处的一阶导数值,返回时存放n个给定结点处的一阶导数值。
    // ddy   : 双精度实型一维数组,长度为n,返回时存放n个给定结点处的二阶导数值。
    //          t     : 双精度实型一维数组,长度为m,存放指定m个插值点的值x(i);要求x(0)<x(1)<....<x(m-1)。
    //          m     : 整形变量,指定插值点的个数, 
    // z     : 双精度实型一维数组,长度为m,返回时存放m个指定插值点处的函数值。   
    //          dz    : 双精度实型一维数组,长度为m,返回时存放m个指定插值点处的一阶导数值。  
    //          ddz   : 双精度实型一维数组,长度为m,返回时存放m个指定插值点处的二阶导数值。
    //【返回值】double型,f(x)从x[0]到x[n-1]的定积分值
    //【用  途】(第一边界条件)三次样条插值
    //【说  明】
    // 
    //***********************************************************************************************/
    double CMathCalculate::Espl1(double *x, double *y, 
      int n, double *dy, 
      double *ddy, double *t,
      int m, double *z, 
      double *dz, double *ddz)
    {
    int i,j;
    double h0,h1,alpha,beta,g,*s;
    s=(double *)malloc(n*sizeof(double));
    s[0]=dy[0]; dy[0]=0.0;
    h0=x[1]-x[0];
    for (j=1;j<=n-2;j++)
    {
    h1=x[j+1]-x[j];
    alpha=h0/(h0+h1);
    beta=(1.0-alpha)*(y[j]-y[j-1])/h0;
    beta=3.0*(beta+alpha*(y[j+1]-y[j])/h1);
    dy[j]=-alpha/(2.0+(1.0-alpha)*dy[j-1]);
    s[j]=(beta-(1.0-alpha)*s[j-1]);
    s[j]=s[j]/(2.0+(1.0-alpha)*dy[j-1]);
    h0=h1;
    }
    for (j=n-2;j>=0;j--)
    dy[j]=dy[j]*dy[j+1]+s[j];
    for (j=0;j<=n-2;j++)
    s[j]=x[j+1]-x[j];
    for (j=0;j<=n-2;j++)
    {
    h1=s[j]*s[j];
    ddy[j]=6.0*(y[j+1]-y[j])/h1-2.0*(2.0*dy[j]+dy[j+1])/s[j];
    }
    h1=s[n-2]*s[n-2];
    ddy[n-1]=6.*(y[n-2]-y[n-1])/h1+2.*(2.*dy[n-1]+dy[n-2])/s[n-2];
    g=0.0;
    for (i=0;i<=n-2;i++)
    {
    h1=0.5*s[i]*(y[i]+y[i+1]);
    h1=h1-s[i]*s[i]*s[i]*(ddy[i]+ddy[i+1])/24.0;
    g=g+h1;
    }
    for (j=0;j<=m-1;j++)
    {
    if (t[j]>=x[n-1]) 
    i=n-2;
    else
    {
    i=0;
    while (t[j]>x[i+1]) i=i+1;
    }
    h1=(x[i+1]-t[j])/s[i];
    h0=h1*h1;
    z[j]=(3.0*h0-2.0*h0*h1)*y[i];
    z[j]=z[j]+s[i]*(h0-h0*h1)*dy[i];
    dz[j]=6.0*(h0-h1)*y[i]/s[i];
    dz[j]=dz[j]+(3.0*h0-2.0*h1)*dy[i];
    ddz[j]=(6.0-12.0*h1)*y[i]/(s[i]*s[i]);
    ddz[j]=ddz[j]+(2.0-6.0*h1)*dy[i]/s[i];
    h1=(t[j]-x[i])/s[i];
    h0=h1*h1;
    z[j]=z[j]+(3.0*h0-2.0*h0*h1)*y[i+1];
    z[j]=z[j]-s[i]*(h0-h0*h1)*dy[i+1];
    dz[j]=dz[j]-6.0*(h0-h1)*y[i+1]/s[i];
    dz[j]=dz[j]+(3.0*h0-2.0*h1)*dy[i+1];
    ddz[j]=ddz[j]+(6.0-12.0*h1)*y[i+1]/(s[i]*s[i]);
    ddz[j]=ddz[j]-(2.0-6.0*h1)*dy[i+1]/s[i];
    }
    free(s);
    return(g);
    }这个方法应该是能解决你的问题!因为这个方法需要知道第一个点和最后一个点的一阶导数值,如果你实在不知道,可以用第二个点的函数值减去第一个点的函数值,再除以第二个结点和第一个结点的x值的差值来作为第一个点的导数值来近似表示,同样用最后一个点的函数值减去倒数第二个点的函数值,再除以最后一个点和倒数第二个点的x插值最为最后一个点的导数值来近似表示,因为“导数”本来就是“微商”嘛!如果你还不明白,请及时回复。