TSpectrum2 deconvolution

Hi,

I am trying to execute a script in order to deconvolve my spectrum in 2 dimensions like in the example 6 : https://root.cern.ch/doc/v608/classTSpectrum2.html. But I got a message that I have a segmentation violation.

I analyse the source code of deconvolution (https://root.cern.ch/doc/v608/TSpectrum2_8cxx_source.html), at the line 1222, it look like the tab working_space has not enough space to stock the data. I solve this problem by increasing the size of the tab working_space at the line 1149-1150 and it works now.

I would like to know if you have the same problem.

Thank.

The size of working_space is defined by ssizex which is an input parameter. What have you done exactly ? can you provide an example ?

I just copy the example 6 and apply to my problem. I have a PSF.txt and a source.txt, I try to execute directly the Deconvolution function but I have a segmentation violation. So I used directly the function source and change the size of the tab working_space.

Archive.zip (6.8 KB)

Thanks for your code but it is not clear to me how it is related to the lines 1222 and 1149 in the orignal code of TSpectrum2.cxx. I thought you had found a bug in that file and wanted us to check/apply a fix you have.

At the line 1149-1150 in the orignal code of TSpectrum2.cxx , the program declares the tab working_space with the size X and Y. I execute the program, and I find a segmentation. I try debug by putting some “printf” in the loop, and at the line 1222 in the orignal code of TSpectrum2.cxx, the segmentation appears here : working_space[i1 - i1min][i2 - i2min + 2 * ssizey ].

The value [i1 - i1min] and [i2 - i2min + 2 * ssizey] exceed the size X and Y.

I hope you understand this.

Yes I understand. But you also said:

My question is: how have you increased the space ? can you show me the fix you did in TSpectrum2.cxx ?

In my code, at the line 210, I increase it, but I don’t know if this fix could have an impact in the calculation.

In the original code, we have :
Double_t **working_space = new Double_t *[ssizex];
for (i = 0; i < ssizex; i++)
working_space[i] = new Double_t[5 * ssizey];

I change it :
Double_t **working_space = new Double_t [2 * ssizex];
for (i = 0; i < 2
ssizex; i++)
working_space[i] = new Double_t[10 * ssizey];

I am not the author of this code so I cannot give you an absolute answer.Can you send me the original small macro you wrote producing the Set Fault ? then I can try it with my ROOT and see I see it also.

Deconvolution2D_simple.C (10.8 KB)

You need the source.txt and PSF_100.txt in my previous answer.

At the line 168-169 you have the 1st method, you execute directly the function, so you need to comment the line 174 to 341.

Or you can execute the source code of function “Deconvolution”, so you need to comment the line 168-169.

The example you sent works for me. I am interested in the one which crashes. Can you provide it ?

This code crashes by using the class Deconvolution.
Deconvolution2D_simple_crash.C (10.8 KB)

Thanks. Yes it crakes for me also.

Thank you for your help.

I hope you can contact the authors of this code TSpectrum2.cxx.

I can’t … he is not available anymore … I try to see what could be wrong… right now I see nothing. May be something wrong in your code ?

I think there is a serious problem with this TSpectrum2::Deconvolution … if you look at the examples they instantiate TSpectrum not TSpectrum2… that’s obviously wrong… So I am afraid this method is not tested.

related post: Problem with TSpectrum2 Deconvolution

Have you got the file TSpectrum2.root directly?

We can also test this example : https://root.cern.ch/doc/v608/classTSpectrum2.html

void Decon2() {
   Int_t i, j;
   Double_t nbinsx = 64;
   Double_t nbinsy = 64;
   Double_t xmin = 0;
   Double_t xmax = (Double_t)nbinsx;
   Double_t ymin = 0;
   Double_t ymax = (Double_t)nbinsy;
   Double_t** source = new Double_t*[nbinsx];
   for (i=0;i<nbinsx;i++)
      source[i]=new Double_t[nbinsy];
   TH2F *decon = new TH2F("decon","Gold deconvolution",nbinsx,xmin,xmax,nbinsy,ymin,ymax);
   TFile *f = new TFile("spectra2/TSpectrum2.root");
   decon=(TH2F*) f->Get("decon2;1");
   Double_t** response = new Double_t*[nbinsx];
   for (i=0;i<nbinsx;i++)
      response[i]=new Double_t[nbinsy];
   TH2F *resp = new TH2F("resp","Response matrix",nbinsx,xmin,xmax,nbinsy,ymin,ymax);
   resp=(TH2F*) f->Get("resp2;1");
   TCanvas *Deconvol = new TCanvas("Deconvolution","Gold deconvolution",10,10,1000,700);
   TSpectrum *s = new TSpectrum();
   for (i = 0; i < nbinsx; i++){
      for (j = 0; j < nbinsy; j++){
         source[i][j] = decon->GetBinContent(i + 1,j + 1);
      }
   }
   for (i = 0; i < nbinsx; i++){
      for (j = 0; j < nbinsy; j++){
         response[i][j] = resp->GetBinContent(i + 1,j + 1);
      }
   }
   s->Deconvolution(source,response,nbinsx,nbinsy,1000,1,1);
   for (i = 0; i < nbinsx; i++){
      for (j = 0; j < nbinsy; j++)
         decon->SetBinContent(i + 1,j + 1, source[i][j]);
   }
   decon->Draw("SURF");
}

No :frowning: … it is lost … anyway the code is wrong … it uses TSpectrum … not TSpectrum2

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