1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 /************************************
17 * AliPMD offline calibration task
19 * Satyajit Jena, IIT Bombay
21 * Fri Feb 12 13:30:19 IST 2010
23 ************************************/
31 #include "AliAnalysisTask.h"
32 #include "AliAnalysisManager.h"
33 #include "AliVEvent.h"
35 #include "AliESDEvent.h"
36 #include "AliAODEvent.h"
37 #include "TObjArray.h"
39 #include "AliPMDOfflineCalibTask.h"
41 ClassImp(AliPMDOfflineCalibTask)
43 //________________________________________________________________________
44 AliPMDOfflineCalibTask::AliPMDOfflineCalibTask(const char *name)
45 : AliAnalysisTaskSE(name),
52 fSelectedTrigger(new TObjArray()),
57 DefineOutput(1, TList::Class());
60 //________________________________________________________________________
61 void AliPMDOfflineCalibTask::UserCreateOutputObjects()
63 fListOfHistos = new TList();
65 fNEvents = new TH1I("hNEvents","Events Statistic", 3, 0, 5);
67 fPmdCalibAdcP = new TH1F("fPmdCalibAdcP", "PMD Calibration ADC fill for PRE", 234795, 0, 234795);
68 fPmdCalibAdcP->GetYaxis()->SetTitle("ADC");
69 fPmdCalibAdcP->GetXaxis()->SetTitle("Cells in SMN-ROW-COL");
71 fPmdCalibAdcC = new TH1F("fPmdCalibAdcC", "PMD Calibration ADC fill for CPV", 234795, 0, 234795);
72 fPmdCalibAdcC->GetYaxis()->SetTitle("ADC");
73 fPmdCalibAdcC->GetXaxis()->SetTitle("Cells in SMN-ROW-COL");
75 fPmdCalibEntP = new TH1F("fPmdCalibEntP", "PMD Calibration Entries fill for PRE", 234795, 0, 234795);
76 fPmdCalibEntP->GetYaxis()->SetTitle("Frequescy");
77 fPmdCalibEntP->GetXaxis()->SetTitle("Cells in SMN-ROW-COL");
79 fPmdCalibEntC = new TH1F("fPmdCalibEntC", "PMD Calibration Entries fill for CPV", 234795, 0, 234795);
80 fPmdCalibEntC->GetYaxis()->SetTitle("Frequescy ");
81 fPmdCalibEntC->GetXaxis()->SetTitle("Cells in SMN-ROW-COL");
83 fListOfHistos->Add(fPmdCalibAdcP);
84 fListOfHistos->Add(fPmdCalibAdcC);
85 fListOfHistos->Add(fPmdCalibEntP);
86 fListOfHistos->Add(fPmdCalibEntC);
87 fListOfHistos->Add(fNEvents);
90 //________________________________________________________________________
91 void AliPMDOfflineCalibTask::UserExec(Option_t *)
94 AliVEvent *event = InputEvent();
96 Printf("ERROR: Could not retrieve event");
102 AliESDEvent* pmdesd = dynamic_cast<AliESDEvent*>(event);
107 Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
112 for(Int_t k = 0; k < numberOfTriggerSelected; k++)
114 const TObjString *const obString = (TObjString*)fSelectedTrigger->At(k);
116 const TString tString = obString->GetString();
117 if(pmdesd->IsTriggerClassFired((const char*)tString))
126 for(Int_t k = 0; k < numberOfTriggerSelected; k++)
128 const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
130 const TString tString = obString->GetString();
131 if(pmdesd->IsTriggerClassFired((const char*)tString))
140 PostData(1, fListOfHistos);
146 Int_t npmdcl = pmdesd->GetNumberOfPmdTracks();
147 if(npmdcl < 1) fNEvents->Fill(4);
151 AliESDPmdTrack *pmdtr = pmdesd->GetPmdTrack(npmdcl);
153 Int_t det = pmdtr->GetDetector();
154 Int_t smn = pmdtr->GetSmn();
155 Float_t adc = pmdtr->GetClusterADC();
156 Float_t sTag = pmdtr->GetClusterSigmaX();
157 Float_t sRowCol = pmdtr->GetClusterSigmaY();
160 Float_t rc = smn*10000 + sRowCol;
161 if(adc > 1200.) continue;
163 if(sTag > 999. && sTag < 1000.)
167 fPmdCalibAdcP->Fill(rc,adc);
168 fPmdCalibEntP->Fill(rc);
172 fPmdCalibAdcC->Fill(rc,adc);
173 fPmdCalibEntC->Fill(rc);
179 PostData(1, fListOfHistos);
182 //________________________________________________________________________
183 void AliPMDOfflineCalibTask::Terminate(Option_t *)
185 fListOfHistos = dynamic_cast<TList*>(GetOutputData(1));
187 fPmdCalibAdcP = dynamic_cast<TH1F*>(fListOfHistos->At(0));
189 printf("ERROR: No ADC File Generated for PMD-PRE ");
193 fPmdCalibAdcC = dynamic_cast<TH1F*>(fListOfHistos->At(1));
195 printf("ERROR: No ADC File Generated for PMD-CPV ");
199 fPmdCalibEntP = dynamic_cast<TH1F*>(fListOfHistos->At(2));
201 printf("ERROR: No Entry File Generated for PMD-PRE ");
202 printf(" No fhXyPRE ");
206 fPmdCalibEntC = dynamic_cast<TH1F*>(fListOfHistos->At(3));
208 printf("ERROR: No Entry File Generated for PMD-CPV ");
212 fNEvents = dynamic_cast<TH1I*>(fListOfHistos->At(4));
214 printf("ERROR: No Entry File Generated for Event Counter ");
219 Info("AliPMDOfflineCalibTask","PMD offline Calibration Successfully finished");