_ROOT Version:6.14/04

Hello everyone,

I try to use the Deconvolution of the TSpectrum class and I get results I cannot understand.
This deconvolution process is based on the gold deconvolution algorithm. You can see explanations in the following website :
Gold deconvolution process.

Deconvolution(double *source,double *response, int size, int numberIterations, int numberRepetitions)
What I understood is that we set the source in argument is the y, the response is the h and after the process, the y is transformed to be x the solution of the system.

I tried to use this to solve a convolution between a gaussian signal (the h) and a random x (that I chose). After deconvolution, I cannot get the x at all. Here is my very simple code :

void Deconv_test() {
/*** Fill the x and h ***/

// x
int N=20;
double *x = new double[N];
x[0]=1.0; x[1]=33.0; x[2]=5.0; x[3]=4.0; x[4]=10.1; x[5]=22.3; x[6]=6.0; x[7]=8.7; x[8]=9.3; x[9]=1.1;
x[10]=0.2; x[11]=0.3; x[12]=1.05; x[13]=5.7; x[14]=1.2; x[15]=0.2; x[16]=3.6; x[17]=7.0; x[18]=1.0; x[19]=0.1;

// elec
int Nr=14.0; 
double sigma = 1.5;
double *elec=new double[Nr];
for(int i = 0 ; i < Nr ; i++) {
	elec[i] = exp(-0.5 * (i - Nr/2) * (i - Nr/2) / sigma / sigma);
        elec[i] = elec[i] / sqrt(2 * M_PI) / sigma;}
double sum=0.0;
for (int i=0 ; i < Nr ; i++) sum += elec[i];
for (int i = 0 ; i < Nr ; i++) elec[i] = elec[i] * 100. / sum;

// h (to be at the same format as the Deconvolution method) OR response !
double *h=new double[N];
for(int i=0;i<N;i++) h[i]=0.0;
for(int i=0;i<N;i++) h[i]=elec[i];

/*** Convolution ***/

double** H = new double*[2*N];  
for (int i=0;i<2*N;i++) H[i]=new double[N];
for (int i = 0; i < 2*N; i++) {
	for (int j = 0; j < N; j++) H[i][j] = 0.0;   }
for (int j = 0; j < N; j++) {
   	for (int i = j; i < j+N; i++) H[i][j] = h[i-j]; }

double *source=new double[2*N]; // or y 
for(int k=0;k<2*N;k++) {
	for(int j=0;j<N;j++) 
		source[k]+=H[k][j]*x[j]; }

/*** Deconvolution ***/

TSpectrum *s = new TSpectrum();
for(int i=0;i<N;i++) cout << x[i] <<" "<< source[i]<<endl;

delete[] source;
delete[] h;
delete[] x;
delete[] elec;
delete[] H;


I printed values of the source after the process and the values of x, they are completely different. There should be something I did not understand in the use of this method.

Thank you in advance for your answer,


I do not know if it helps but the documentation is here:

and the source code there:

Thank you for your answer Couet.
I have ever seen these documentations but it did not help to find the problem unfortunately.

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.