TMarLumi.C

Go to the documentation of this file.

// Class TMarLumi
//
// Author     : F. Cassol Brunner
// Created    : 02/2002
// Last update: $Date: 2005/06/24 08:35:29 $
//          by: $Author: gbrandt $
// Comment: Contains definition of lumi class

#ifndef __TMARLUMI
#include "Marana/TMarLumi.h"
#endif

#define MARLUMIDEBUG 0



ClassImp(TMarLumi)

TMarLumi::TMarLumi(string filename, string FileBadRun="not_given")
{


 // Init Lumi
 RunNumber=0;
 RunNumberOld=0;
 LumiPerRun=0;
 TotalLumi=0;
 LumiPerRunOld=0;
 TotalLumiOld=0;
 detStatus= NULL; 
  runList =NULL;  
//...change  TMarLumi::GetPolar() if you change the bin of these histograms!!!!
 TH1F* hpolarR = new TH1F("polarR","polarisation Right",33, -82.5, 82.5);
 TH1F* hpolarL = new TH1F("polarL","polarisation Left",33, -82.5, 82.5);
 TH1F* hpolar  = new TH1F("polar","polarisation All",33, -82.5, 82.5);
 
 polar_integral = NULL;

 float nLzero=0;
 float nRzero=0;



 rpolar     = new TRandom();

 try{

   Int_t nrun;
   Float_t lumi;
   Float_t polar=0;

   // Save current ROOT directory so histos below get added to correct dir
   TDirectory *savedir = gDirectory;
   TFile* runfile=new TFile(filename.c_str());
   savedir->cd();

   if (!runfile){
       Error("TMarLumi","Cannot open lumi file '%s'",runfile->GetName());
       throw string("exit");
   }
   runfile->GetObject("H1RunList",runList);
   if (!runList) {
        Error("TMarLumi","Cannot get H1RunList from '%s'",runfile->GetName());
        throw  string("exit");
   }
   runfile->GetObject("MyDetectorStatus",detStatus);
   if (!detStatus) {
        Error("TMarLumi","Cannot get H1DetectorStatus from '%s'",runfile->GetName());
        throw  string("exit");
   }
   else { detStatus->Print(); }

   runList->Sort();

   for(int n=0;n<runList->GetNrOfRuns();n++){
       H1RunLumi* runLumi=(H1RunLumi*)runList->At(n);

       nrun=runLumi->GetRunNumber();
       lumi=runLumi->GetCorrectedLumi();
       polar=runLumi->GetPolarisation();

       if(MARLUMIDEBUG)cout << "- good run "<< nrun << " lumi "<< lumi << " polarization "<< polar<<endl;

       if(IsDataR(nrun)==kTRUE){
         if(polar!=0)hpolarR->Fill(polar,lumi);
          else nRzero+=lumi;
       }
       else if(IsDataL(nrun)==kTRUE){
          if(polar!=0)hpolarL->Fill(polar,lumi);
          else nLzero+=lumi;
       }
   }
   runfile->Close();
   delete  runfile;

   printf("- Found %d good runs in Lumi file \n",runList->GetNrOfRuns());



  } //... try
  catch(string error)
 {
     cout << error <<endl;;
     exit(0); 
  
 }
 catch(...)
 {
    throw string("TMarLumi::TMarLumi = Undefined error in reading parameters_file! - stop");
    exit(0); 
 }     


//...scale the left and right histogram adding the run without pol data
  
  if(nRzero>0&&hpolarR->Integral()>0){
    float factor=(hpolarR->Integral()+nRzero)/hpolarR->Integral();
    hpolarR->Scale(factor);
  }  
  if(nLzero>0&&hpolarL->Integral()>0){
    float factor=(hpolarL->Integral()+nLzero)/hpolarL->Integral();
    hpolarL->Scale(factor);
  }

  hpolar->Add(hpolarR);
  hpolar->Add(hpolarL);
  
   
//... vector with  profile of the polarisation integral (normalized to 1)
  
  if(hpolar->GetEntries()>0){
    hpolar->ComputeIntegral();
    polar_integral = new double[33];
    polar_integral=hpolar->GetIntegral();
  }


}




// Standard destructor
TMarLumi::~TMarLumi()
{

};

  
Bool_t TMarLumi::Update(Int_t nrun)
{
//     Get the lumi of the run number

 RunNumberOld=RunNumber;
 LumiPerRunOld=LumiPerRun;
 TotalLumiOld=TotalLumi;

//...search in the list of good lumi runs
 if(runList->IsSelected(nrun)){

     LumiPerRun=runList->GetIntLumi(nrun,nrun);

 }
 else
 {
//...else search in the list of bad lumi runs

  //cout <<"--> TMarLumi::Update = Run "<< nrun <<" is not in the RunList"<<endl;

   return(kFALSE);
 }


 TotalLumi+=LumiPerRun;
 RunNumber=nrun;

 return(kTRUE);  
 
 
}

Float_t TMarLumi::GetPolar(){

  if(polar_integral==NULL) return 0;

//... get a uniform random number between 0 and 1
  Double_t u=rpolar->Uniform();
  int nbin=TMath::BinarySearch(33,polar_integral,u);
  
  Float_t pola=(-80.+nbin*5)/100; 

  return pola;
};

Bool_t TMarLumi::IsDataL( Int_t& irun)
{
  if((irun > 376562  && irun <= 386696 )     || //0304L1
     (irun >= 402992 && irun  <=411192))        //05L1
    {
      
      return kTRUE;
    }
  return kFALSE;
};
Bool_t TMarLumi::IsDataR( Int_t& irun)
{
  if( (irun >= 357072 && irun <= 376562) ||    //0304R1
      (irun >= 387942 && irun <= 392213) ||    //0304R2
      (irun >= 399101 && irun <= 402634)  ||   //05R1
      (irun >  411192))                        //05R2
  {
     return kTRUE;
  }
  return kFALSE;
}

---