TRD/qaRec/AliTRDrecoTask.cxx

   1 #include "TClass.h"
   2 #include "TMethod.h"
   3 #include "TMethodCall.h"
   4 #include "TMethodArg.h"
   5 #include "TFile.h"
   6 #include "TList.h"
   7 #include "TH1.h"
   8 #include "TF1.h"
   9 #include "TObjArray.h"
  10 #include "TDirectory.h"
  11 #include "TTreeStream.h"
  12
  13 #include "AliLog.h"
  14 #include "AliAnalysisTask.h"
  15
  16 #include "AliTRDrecoTask.h"
  17
  18 ClassImp(AliTRDrecoTask)
  19 FILE* AliTRDrecoTask::fgFile = 0x0;
  20 //_______________________________________________________
  21 AliTRDrecoTask::AliTRDrecoTask(const char *name, const char *title)
  22   : AliAnalysisTask(name, title)
  23   ,fNRefFigures(0)
  24   ,fDebugLevel(0)
  25   ,fPlotFuncList(0x0)
  26   ,fContainer(0x0)
  27   ,fTracks(0x0)
  28   ,fTrack(0x0)
  29   ,fMC(0x0)
  30   ,fESD(0x0)
  31   ,fDebugStream(0x0)
  32 {
  33   DefineInput(0, TObjArray::Class());
  34   DefineOutput(0, TObjArray::Class());
  35 }
  36
  37 //_______________________________________________________
  38 AliTRDrecoTask::~AliTRDrecoTask()
  39 {
  40   printf(" %s (%s)\n", GetName(), GetTitle());
  41   if(fDebugStream){
  42     delete fDebugStream;
  43     fDebugStream = 0x0;
  44   }
  45
  46   if(fPlotFuncList){
  47     fPlotFuncList->Delete();
  48     delete fPlotFuncList;
  49     fPlotFuncList = 0x0;
  50   }
  51
  52   if(fContainer){
  53     if(fContainer->IsOwner()) fContainer->Delete();
  54     delete fContainer;
  55     fContainer = 0x0;
  56   }
  57
  58   if(fgFile){
  59     fflush(fgFile);
  60     fclose(fgFile);
  61   }
  62 }
  63
  64 //_______________________________________________________
  65 void AliTRDrecoTask::ConnectInputData(Option_t *)
  66 {
  67   //
  68   // Connect input data
  69   //
  70
  71   fTracks = dynamic_cast<TObjArray *>(GetInputData(0));
  72 }
  73
  74 //_______________________________________________________
  75 void AliTRDrecoTask::Exec(Option_t *)
  76 {
  77   if(!fPlotFuncList){
  78     AliWarning("No functor list defined for the reference plots");
  79     return;
  80   }
  81   if(!fTracks) return;
  82   if(!fTracks->GetEntriesFast()) return;
  83
  84   AliTRDtrackInfo *trackInfo = 0x0;
  85   TIter plotIter(fPlotFuncList);
  86   TObjArrayIter trackIter(fTracks);
  87   while((trackInfo = dynamic_cast<AliTRDtrackInfo*>(trackIter()))){
  88     fTrack = trackInfo->GetTrack();
  89     fMC    = trackInfo->GetMCinfo();
  90     fESD   = trackInfo->GetESDinfo();
  91
  92     TMethodCall *plot = 0x0;
  93     plotIter.Reset();
  94     while((plot=dynamic_cast<TMethodCall*>(plotIter()))){
  95       plot->Execute(this);
  96     }
  97   }
  98   PostData(0, fContainer);
  99 }
 100
 101 //_______________________________________________________
 102 Bool_t AliTRDrecoTask::GetRefFigure(Int_t /*ifig*/)
 103 {
 104   AliWarning("Retrieving reference figures not implemented.");
 105   return kFALSE;
 106 }
 107
 108 //_______________________________________________________
 109 Bool_t AliTRDrecoTask::PutTrendValue(Char_t *name, Double_t val, Double_t err)
 110 {
 111   if(!fgFile){
 112     fgFile = fopen("TRD.Performance.txt", "wt");
 113   }
 114   fprintf(fgFile, "%s_%s %f %f\n", GetName(), name, val, err);
 115   return kTRUE;
 116 }
 117
 118 //_______________________________________________________
 119 void AliTRDrecoTask::InitFunctorList()
 120 {
 121   TClass *c = this->IsA();
 122
 123   TMethod *m = 0x0;
 124   TIter methIter(c->GetListOfMethods());
 125   while((m=dynamic_cast<TMethod*>(methIter()))){
 126     TString name(m->GetName());
 127     if(!name.BeginsWith("Plot")) continue;
 128     if(!fPlotFuncList) fPlotFuncList = new TList();
 129     fPlotFuncList->AddLast(new TMethodCall(c, (const char*)name, ""));
 130   }
 131 }
 132
 133 //_______________________________________________________
 134 Bool_t AliTRDrecoTask::Load(const Char_t *filename)
 135 {
 136   if(!TFile::Open(filename)){
 137     AliWarning(Form("Couldn't open file %s.", filename));
 138     return kFALSE;
 139   }
 140   TObjArray *o = 0x0;
 141   if(!(o = (TObjArray*)gFile->Get(GetName()))){
 142     AliWarning("Missing histogram container.");
 143     return kFALSE;
 144   }
 145   fContainer = (TObjArray*)o->Clone(GetName());
 146   gFile->Close();
 147   return kTRUE;
 148 }
 149
 150 //________________________________________________________
 151 Bool_t AliTRDrecoTask::Save(TObjArray *results){
 152   //
 153   // Store the output graphs in a ROOT file
 154   // Input TObject array will not be written as Key to the file,
 155   // only content itself
 156   //
 157
 158   TDirectory *cwd = gDirectory;
 159   if(!TFile::Open(Form("TRD.Result%s.root", GetName()), "RECREATE")) return kFALSE;
 160
 161   TIterator *iter = results->MakeIterator();
 162   TObject *inObject = 0x0, *outObject = 0x0;
 163   while((inObject = iter->Next())){
 164     outObject = inObject->Clone();
 165     outObject->Write(0x0, TObject::kSingleKey);
 166   }
 167   delete iter;
 168   gFile->Close(); delete gFile;
 169   cwd->cd();
 170   return kTRUE;
 171 }
 172
 173 //_______________________________________________________
 174 Bool_t AliTRDrecoTask::PostProcess()
 175 {
 176   AliWarning("Post processing of reference histograms not implemented.");
 177   return kFALSE;
 178 }
 179
 180 //_______________________________________________________
 181 void AliTRDrecoTask::SetDebugLevel(Int_t level)
 182 {
 183   fDebugLevel = level;
 184   if(fDebugLevel>=1){
 185     TDirectory *savedir = gDirectory;
 186     fDebugStream = new TTreeSRedirector(Form("TRD.DBG%s.root", GetName()));
 187     savedir->cd();
 188   }
 189 }
 190
 191 //________________________________________________________
 192 void AliTRDrecoTask::Adjust(TF1 *f, TH1 *h)
 193 {
 194 // Helper function to avoid duplication of code
 195 // Make first guesses on the fit parameters
 196
 197   // find the intial parameters of the fit !! (thanks George)
 198   Int_t nbinsy = Int_t(.5*h->GetNbinsX());
 199   Double_t sum = 0.;
 200   for(Int_t jbin=nbinsy-4; jbin<=nbinsy+4; jbin++) sum+=h->GetBinContent(jbin); sum/=9.;
 201   f->SetParLimits(0, 0., 3.*sum);
 202   f->SetParameter(0, .9*sum);
 203
 204   f->SetParLimits(1, -.2, .2);
 205   f->SetParameter(1, -0.1);
 206
 207   f->SetParLimits(2, 0., 4.e-1);
 208   f->SetParameter(2, 2.e-2);
 209   if(f->GetNpar() <= 4) return;
 210
 211   f->SetParLimits(3, 0., sum);
 212   f->SetParameter(3, .1*sum);
 213
 214   f->SetParLimits(4, -.3, .3);
 215   f->SetParameter(4, 0.);
 216
 217   f->SetParLimits(5, 0., 1.e2);
 218   f->SetParameter(5, 2.e-1);
 219 }