Help on Method/s to create Tree from multi-fold gamma-ray coincidence data

Wile_E_Coyote · August 11, 2018, 7:08am

I think you can easily reduce the amount of the required RAM by (up to) a factor 6, if you explicitly sort your gammas:

#include <algorithm>
// ...
      for(int j = 0; j < (No_Clover - 2); j++) {
        for(int k = (j + 1); k < (No_Clover - 1); k++) { // (k > j)
          for(int l = (k + 1); l < No_Clover; l++) { // (l > k)
            // ...
                        // Double_t x[3] = {(Double_T_t)Eclab[j], (Double_T_t)Eclab[k], (Double_T_t)Eclab[l]};
                        // Double_t x[3] = {(Eclab[j] + 1.e-6), (Eclab[k] + 1.e-6), (Eclab[l] + 1.e-6)};
                        Double_t x[3] = {(Eclab[j] + 0.5), (Eclab[k] + 0.5), (Eclab[l] + 0.5)};
                        std::sort(x, x + 3); // sort in an increasing order
                        if (his3D->GetNbins() <= MaxNbins) his3D->Fill(x);
                        // ...

When you do not sort your gammas, each triple of energies appears 6 times in your 3D histogram so, there are 6 “peaks” in different areas of the 3D histogram (where 6 is the number of permutations, of course).

If you sort them (in an increasing order), you will get 1 “peak” only (always energy_gamma_1 <= energy_gamma_2 <= energy_gamma_3).

This makes your 3D histogram “totally unsymmetric” but, I think you can take it into account when projecting / analysing it later (well, I’m afraid that this idea requires some detailed thoughtful examination so, you should discuss it with your colleagues).

ajaydeo · August 11, 2018, 11:49am

I am open to make thorough checks using the procedure that you have mentioned. However, I don’t to how to analyze such “unsymmetric” cube. Can you please tell me the possible analysis method/s?

If I have to make such a “totally unsymmetric” cube, then probably, TH3F class will also work. Then, I can make use of TSpectrum3 for subtracting 3D-background from the 3D-histogram, which I guess not possible in THnSparse.

BTW, I am able to successfully generate 3D THnSparse with lower dimensions (1000) with full data set.

ajaydeo · August 11, 2018, 3:47pm

Is this same as scaling or using weightage while filling the histogram?
i.e. his3D->Fill(x,0.1666666);

Wile_E_Coyote · August 12, 2018, 9:08am

Maybe this works:

#include "TH1.h"
#include "TH2.h"
#include "THnSparse.h"
#include "TAxis.h"
#include "TString.h"

// Usage example: TH1D *h1 = Project1D(his3D, 449.5, 452.5, 668.5, 671.5);
// Note: in principle, you should set low1 <= up1 and low2 <= up2
TH1D *Project1D(THnSparse *h, // the ("totally unsymmetric") "cubic" histogram
                Double_t low1, Double_t up1, // the first "window"
                Double_t low2, Double_t up2, // the second "window"
                Option_t *option = "") { // THnSparse::Projection option
  if (!h) return 0; // just a precaution
  TAxis *a0 = h->GetAxis(0);
  TAxis *a1 = h->GetAxis(1);
  TAxis *a2 = h->GetAxis(2);
  a0->SetRange(a0->FindFixBin( low1 ), a0->FindFixBin( up1 ));
  a1->SetRange(a1->FindFixBin( low2 ), a1->FindFixBin( up2 ));
  a2->SetRange(-1, -1);
  TH1D *h1 = h->Projection(2, option);
  h1->SetName(TString::Format("%s_Project1D", h->GetName()));
  a0->SetRange(a0->FindFixBin( low2 ), a0->FindFixBin( up2 ));
  a1->SetRange(a1->FindFixBin( low1 ), a1->FindFixBin( up1 ));
  TH1D *hi = h->Projection(2, option);
  h1->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(a0->FindFixBin( low1 ), a0->FindFixBin( up1 ));
  a1->SetRange(-1, -1);
  a2->SetRange(a2->FindFixBin( low2 ), a2->FindFixBin( up2 ));
  hi = h->Projection(1, option);
  h1->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(a0->FindFixBin( low2 ), a0->FindFixBin( up2 ));
  a2->SetRange(a2->FindFixBin( low1 ), a2->FindFixBin( up1 ));
  hi = h->Projection(1, option);
  h1->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(-1, -1);
  a1->SetRange(a1->FindFixBin( low1 ), a1->FindFixBin( up1 ));
  a2->SetRange(a2->FindFixBin( low2 ), a2->FindFixBin( up2 ));
  hi = h->Projection(0, option);
  h1->Add(hi);
  delete hi; // no longer needed
  a1->SetRange(a1->FindFixBin( low2 ), a1->FindFixBin( up2 ));
  a2->SetRange(a2->FindFixBin( low1 ), a2->FindFixBin( up1 ));
  hi = h->Projection(0, option);
  h1->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(-1, -1);
  a1->SetRange(-1, -1);
  a2->SetRange(-1, -1);
  return h1;
}

// Usage example: TH2D *h2 = Project2D(his3D, 449.5, 452.5);
// Note: in principle, you should set low <= up
TH2D *Project2D(THnSparse *h, // the ("totally unsymmetric") "cubic" histogram
                Double_t low, Double_t up, // the "window"
                Option_t *option = "") { // THnSparse::Projection option
  if (!h) return 0; // just a precaution
  TAxis *a0 = h->GetAxis(0);
  TAxis *a1 = h->GetAxis(1);
  TAxis *a2 = h->GetAxis(2);
  a0->SetRange(a0->FindFixBin( low ), a0->FindFixBin( up ));
  a1->SetRange(-1, -1);
  a2->SetRange(-1, -1);
  TH2D *h2 = h->Projection(1, 2, option);
  h2->SetName(TString::Format("%s_Project2D", h->GetName()));
  TH2D *hi = h->Projection(2, 1, option);
  h2->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(-1, -1);
  a1->SetRange(a1->FindFixBin( low ), a1->FindFixBin( up ));
  a2->SetRange(-1, -1);
  hi = h->Projection(0, 2, option);
  h2->Add(hi);
  delete hi; // no longer needed
  hi = h->Projection(2, 0, option);
  h2->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(-1, -1);
  a1->SetRange(-1, -1);
  a2->SetRange(a2->FindFixBin( low ), a2->FindFixBin( up ));
  hi = h->Projection(0, 1, option);
  h2->Add(hi);
  delete hi; // no longer needed
  hi = h->Projection(1, 0, option);
  h2->Add(hi);
  delete hi; // no longer needed
  a0->SetRange(-1, -1);
  a1->SetRange(-1, -1);
  a2->SetRange(-1, -1);
  return h2;
}

ajaydeo · August 13, 2018, 11:02am

Please give me some time to test this. Will get back to you with the outcome after testing.

In the meantime, can you please answer my last question?

You have mentioned in the post 61 of this thread that we can reduce the RAM usage by a factor of 6 by creating an “asymmetric” 3D-histogram. Is it the same as using

his3D->Fill(x,0.1666666); // 1/6 = 0.16666666

instead of

his3D->Fill(x);

Wile_E_Coyote · August 13, 2018, 11:04am

Just try it.

ajaydeo · August 14, 2018, 10:04am

I tried it, and it failed to work

Now, I am testing the last method. Will get back to you with results.

BTW, do you know whether @pcanal is working on/has time to resolve the issue?

ajaydeo · August 15, 2018, 12:28pm

I have tested the “unsymmetric” cubic histogram and the corresponding projection method that you have suggested. The preliminary observations indicate that everything is working perfectly as expected. I will do more thorough check in the weekend.

As I had mentioned earlier, I have to create several projection1D using different combination of gamma-ray pairs. For this, I am using the macro “Sparse_project.C” (attached). I use the following sequence to generate two 1D-Projections:

root [0] TFile *f = TFile::Open("sparse3D_asym_all_4k.root");
root [1] .ls
TFile**		sparse3D_asym_all_4k.root	
 TFile*		sparse3D_asym_all_4k.root	
  KEY: THnSparseT<TArrayF>	sparse3D;1	sparse3D
root [2] .L Sparse_project.C 
root [3] Project1D(sparse3D,449.5,453.5,668.5,672.5);
root [4] .ls
TFile**		sparse3D_asym_all_4k.root	
 TFile*		sparse3D_asym_all_4k.root	
  OBJ: TH1D	sparse3D_Project1D	sparse3D projection EclabZ : 0 at: 0xec6d440
  KEY: THnSparseT<TArrayF>	sparse3D;1	sparse3D
root [5] sparse3D_Project1D->Draw();
Info in <TCanvas::MakeDefCanvas>:  created default TCanvas with name c1
root [6] sparse3D_Project1D->SetName("his_451_670");
root [7] his_451_670->Draw();
root [8] Project1D(sparse3D,289.5,293.5,698.5,702.5);
root [9] .ls
TFile**		sparse3D_asym_all_4k.root	
 TFile*		sparse3D_asym_all_4k.root	
  OBJ: TH1D	his_451_670	sparse3D projection EclabZ : 0 at: 0xec6d440
  OBJ: TH1D	sparse3D_Project1D	sparse3D projection EclabZ : 0 at: 0xf129740
  KEY: THnSparseT<TArrayF>	sparse3D;1	sparse3D
root [10] sparse3D_Project1D->Draw();
root [11] sparse3D_Project1D->SetName("his_291_700");
root [12] his_291_700->Draw();
root [13] his_451_670->Draw();
root [14]

However, there is a problem. Both the spectra are identical. As you can see, I am using entirely different sets of gamma rays to generate two 1D-histograms. Am I doing anything wrong?

Also, I feel that it would be ideal to mention only gamma-ray energies while projecting, along with a variable “width parameter” which will internally set the “window” for the gamma rays.

Sparse_project.C (3.4 KB)

gen_Sparse3D.cxx (6.9 KB)

Wile_E_Coyote · August 15, 2018, 1:43pm

Try (first you need to root [0] .L Sparse_project.C, of course):

{
  TFile *f = TFile::Open("sparse3D_asym_all_4k.root");
  THnSparse *h; f->GetObject("sparse3D", h);
  Double_t g1 = 451.5;
  Double_t g2 = 670.5;
  Double_t w = 2.0;
  TH1D *h1 = Project1D(h, g1 - w, g1 + w, g2 - w, g2 + w);
  h1->SetName(TString::Format("his_%d_%d", Int_t(g1), Int_t(g2)));
  h1->Print();
  g1 = 291.5;
  g2 = 700.5;
  TH1D *h2 = Project1D(h, g1 - w, g1 + w, g2 - w, g2 + w);
  h2->SetName(TString::Format("his_%d_%d", Int_t(g1), Int_t(g2)));
  h2->Print();
  TCanvas *c = new TCanvas("c", "c");
  c->Divide(1, 2);
  c->cd(1);
  h1->Draw();
  c->cd(2);
  h2->Draw();
  c->cd(0);
  // do we get exactly the same output as above?
  h1->Print();
  h2->Print();
}

ajaydeo · August 16, 2018, 10:06am

It worked, thank you. The two 1D-histograms are no longer identical.

Please check the attached modified version of the macro which I used for 1D- and 2D- projections of THnSparse.

We know that we have TSpectrum2 and TSpectrum3 classes available for background generation of TH2 and TH3 histograms, respectively.

Are there any equivalent “background generation” methods/procedures available for THnSparse Class?

Sparse_project.C (3.7 KB)

Wile_E_Coyote · August 16, 2018, 10:47am

You don’t need to care if “eg_low <= eg_hi”.
Use “SetName” once per histogram (just modify the first occurrence).

system · August 30, 2018, 10:47am

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