DoGammaFilter.C

Go to the documentation of this file.

TH1F *DoDoppler(Axis_t min_e = 0,  Axis_t max_e = 4096)
{
        // Doppler shift
        Float_t beta = 0.02; 
        // Float_t beta = 0; 
        
        // opening angle 
        Float_t dtheta = 10.0 * TMath::Pi() / 180; // euroball
        // Float_t dtheta = 1 * TMath::Pi() / 180;      // AGATA 
        
        // intrinseque resolution
        // linear curve so two points are needed
        // e1, e2 - energies in keV
        // r1, r2 - half-max width intrinseque
        Float_t e1 = 100.0, r1 = 1.0, e2 = 1300.0, r2 = 2.5;
        
        //now determine the linear function 
        TF1 *f_de_int = new TF1("resolution_int","pol1",min_e,max_e);
        
        f_de_int->SetParameters(r1-e1*(r2-r1)/(e2-e1),(r2-r1)/(e2-e1));
        //f_de_int->Draw();
        
        Int_t bin = Int_t ( max_e -  min_e );
        TH1F *h = new TH1F("GammaFilter_Doppler","Doppler",bin,min_e,max_e);
        
        Float_t de_int, de_theta;
        for (int i = 0; i < bin; i++ ) {
                de_int   = f_de_int->Eval(min_e+i);
                de_theta = i * beta * dtheta ;

                // printf("%d %f %f \n",i,de_int,de_theta);
                h->SetBinContent(i+1,TMath::Sqrt( de_int*de_int + de_theta*de_theta ) / 2.3548 ); // 2.3548 for sigma
        }
        
        h->SetLineColor(2);     
        h->SetBinContent(0,0); h->SetBinContent(bin+1,0); 
        
        h->Draw();
        
        return h;
}

/*
TH1F *DoDoppler(Axis_t min_e = 0,  Axis_t max_e = 4096)
{
        // linear curve so two points are needed
        // e1, e2 - energies in keV
        // r1, r2 - half-max width ..
        Float_t e1 = 100.0, r1 = 1.0, e2 = 1300.0, r2 = 2.5;
        
        // change r1 and r2 into sigma1, sigma2 .. 
        r1 =  r1 / 2.3548; r2 = r2 / 2.3548;
        
        //now determine the linear function 
        TF1 *f = new TF1("resolution","pol1",min_e,max_e);
        
        f->SetParameters(r1-e1*(r2-r1)/(e2-e1),(r2-r1)/(e2-e1));
        //f->Draw();
        
        Int_t bin = Int_t ( max_e -  min_e );
        TH1F *h = new TH1F("GammaFilter_Doppler","Doppler",bin,min_e,max_e);
        
        for (int i = 0; i < bin; i++ ) {
                h->SetBinContent(i+1,f->Eval(min_e+i));
        }
        
        h->SetLineColor(2);     
        h->SetBinContent(0,0); h->SetBinContent(bin+1,0); 
        
        h->Draw();
        
        return h;
}
*/

TH1F *DoCompton(Float_t pt = 0.5, Axis_t min_e = 0,  Axis_t max_e = 4096)
{
        Int_t cut1 = 100, cut2 = 300; Float_t a, b;
        
        a = (1.0-pt)/(cut1-cut2); b = 1.0 - cut1 * a;
        
        Int_t bin = Int_t ( max_e -  min_e );
        TH1F *h = new TH1F("GammaFilter_Compton","Compton",bin,min_e,max_e);
        
        for (int i = 0; i < bin; i++ ) { 
                if ( i < cut1 ) {
                        h->SetBinContent(i+1,1); 
                }
                if ( i >= cut1 && i <= cut2 ) {
                        h->SetBinContent(i+1,a*i+b); 
                }       
                if ( i > cut2 ) {                       
                        h->SetBinContent(i+1,pt); 
                }
        }
        
        h->SetLineColor(2);     
        h->SetBinContent(0,0); h->SetBinContent(bin+1,0); 
        
        h->Draw();
        
        return h;
}