Hilbert(希尔伯特)变换的算法
关键词：
希尔伯特
/******************************************/
/* main.c
/******************************************/
#include
#include
#define N 500
#define sq(X) ((X)*(X))
int main()
{
int i;
double x;
double delta[N];
double kappa[N];
double y;
double xmin = -10.;
double xmax = 10.;
double cd = -1.;
double w = 2.;
double h = (xmax - xmin) / N;
for (i = 0; i < N; i++)
{
x = 2. * (xmin + i * h - cd) / w;
if (x < -1.)
delta[i1] = 0.;
else if (x < 1.)
delta[i1] = sqrt(1. - sq(x));
else
delta[i1] = 0.;
}
hilbert(n, delta, kappa);
for (i = 0; i < N; i++)
{
x = 2. * (xmin + i * h - cd) / w;
if (x < -1.)
y = x + sqrt(sq(x) - 1.);
else if (x < 1.)
y = x;
else
y = x - sqrt(sq(x) - 1.);
(void) printf("%.3f %.3f %.3f %.3f\n", xmin + i * h, delta[i], k
ppa[i], y);
}
return 0;
}
/******************************************/
/* hilbert.c
/******************************************/
void hilbert(int, double[], double[]);
/******************************************/
/* hilbert.c
/******************************************/
#include
#include "hilbert.h"
#define PI 3.14159265
void hilbert(int n, double delta[], double kappa[])
{
double d1, d2, d3, d4;
int i1, i2;
for (i1 = 0; i1 < n; i1++)
{
kappa[i1] = 0.;
for (i2 = 1; i2 < n; i2++)
{
d1 = (i1+i2< 0.;
d2 = (i1-i2>=0)? delta[i1-i2]: 0.;
d3 = (i1+i2+1< 0.;
d4 = (i1-i2-1>=0)? delta[i1-i2-1]: 0.;
kappa[i1] -= 0.5 * (d1-d2) / i2 + 0.5 * (d3 - d4) / (i2
1);
}
kappa[i1] /= PI;
}
}
/***********************************************
/* 说明
/***********************************************
The Hilbert transform
This package performs a numerical Hilbert transform.
Download Compressed tar-file
or
hilbert.c and hilbert.h
Header file
The program must include the header file
#include "hilbert.h"
/***********************************************
/* 说明
kappa[i1] -= 0.5 * (d1-d2) / i2 + 0.5 * (d3 - d4) / (i2
1);
}
kappa[i1] /= PI;
}
}
/***********************************************
/* 说明
/***********************************************
The Hilbert transform
This package performs a numerical Hilbert transform.
Download Compressed tar-file
or
hilbert.c and hilbert.h
Header file
The program must include the header file
#include "hilbert.h"
and then call hilbert() to perform the transform.
Usage
If delta and kappa are arrays of n doubles, both arrays are allocated by the mai
program.
The input data are in delta and the call hilbert(n, delta, kappa) will return th
Hilbert transform of delta in kappa.
n and delta are not modified.
Compilation
Compile the package using a C-compiler.
If you use a makefile and this makefile looks like this
a.out: a.o b.o c.o
cc