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 **************************************************************************/
18 // Objects of this class contain basis for absolute calibrations
26 #include "AliEMCALCalibAbs.h"
30 ClassImp(AliEMCALCalibAbs)
32 //____________________________________________________________________________
33 AliEMCALCalibAbs::AliEMCALCalibAbs(const int nSM) :
37 //Default constructor.
38 for (int i=0; i<fNSuperModule; i++) {
39 fSuperModuleData.Add(new AliEMCALSuperModuleCalibAbs(i));
41 fSuperModuleData.Compress(); // compress the TObjArray
42 fSuperModuleData.SetOwner(kTRUE);
45 //____________________________________________________________________________
46 void AliEMCALCalibAbs::ReadTextCalibAbsInfo(Int_t nSM, const TString &txtFileName,
49 //Read data from txt file. ; coordinates given on SuperModule basis
51 std::ifstream inputFile(txtFileName.Data());
53 printf("AliEMCALCalibAbs::ReadCalibAbsInfo - Cannot open the APD info file %s\n", txtFileName.Data());
59 Int_t iSM = 0; // SuperModule index
63 // list of values to be read
64 // first: overall values for the whole SuperModule
65 Int_t iCalibMethod = 0;
67 Float_t absoluteCalib = 0;
68 // third: info for each tower
69 Float_t relativeCalib = 0; // (ADC>GeV relative gain/conversion), value around 1
72 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
74 for (Int_t i = 0; i < fNSuperModule; i++) {
75 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
77 printf("AliEMCALCalibAbs::ReadCalibAbsInfo - Error while reading input file; likely EOF..\n");
81 t->SetSuperModuleNum(iSM);
83 // first: overall values for the whole SuperModule
84 inputFile >> iCalibMethod >> iCalibPass >> absoluteCalib;
85 t->SetCalibMethod(iCalibMethod);
86 t->SetCalibPass(iCalibPass);
87 t->SetAbsoluteCalib(absoluteCalib);
89 // third: info for each tower
90 for (Int_t j=0; j<nAPDPerSM; j++) {
91 inputFile >> iCol >> iRow >> relativeCalib;
93 // check that input values are not out bounds
94 if (iCol<0 || iCol>(AliEMCALGeoParams::fgkEMCALCols-1) ||
95 iRow<0 || iRow>(AliEMCALGeoParams::fgkEMCALRows-1) ) {
96 printf("AliEMCALCalibAbs::ReadCalibAbsInfo - Error while reading input file; j %d iCol %d iRow %d\n", j, iCol, iRow);
100 // assume that this info is already swapped and done for this basis?
102 // C side, oriented differently than A side: swap is requested
103 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
104 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
107 t->SetRelativeCalib(iCol, iRow, relativeCalib);
117 //____________________________________________________________________________
118 void AliEMCALCalibAbs::WriteTextCalibAbsInfo(const TString &txtFileName,
121 // write data to txt file. ; coordinates given on SuperModule basis
123 std::ofstream outputFile(txtFileName.Data());
125 printf("AliEMCALCalibAbs::WriteCalibAbsInfo - Cannot open the APD output file %s\n", txtFileName.Data());
132 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
133 Float_t relativeCalib = 0;
134 for (Int_t i = 0; i < fNSuperModule; i++) {
135 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
137 // first: overall values for the whole SuperModule
138 outputFile << t->GetSuperModuleNum() << endl;
139 outputFile << t->GetCalibMethod() << " "
140 << t->GetCalibPass() << " "
141 << t->GetAbsoluteCalib() << endl;
143 // third: info for each tower
144 for (Int_t j=0; j<nAPDPerSM; j++) {
145 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
146 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
148 relativeCalib = t->GetRelativeCalib(iCol, iRow);
151 // C side, oriented differently than A side: swap is requested
152 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
153 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
156 outputFile << iCol << " " << iRow
157 << " " << relativeCalib << endl;
167 //____________________________________________________________________________
168 void AliEMCALCalibAbs::ReadRootCalibAbsInfo(const TString &rootFileName,
171 //Read data from root file. ; coordinates given on SuperModule basis
172 TFile inputFile(rootFileName, "read");
174 TTree *tree = (TTree*) inputFile.Get("tree");
176 ReadTreeCalibAbsInfo(tree, swapSides);
183 //____________________________________________________________________________
184 void AliEMCALCalibAbs::ReadTreeCalibAbsInfo(TTree *tree,
187 // how many SuperModule's worth of info do we have?
188 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
189 fNSuperModule = tree->GetEntries();
191 Int_t iSM = 0; // SuperModule index
192 // list of values to be read
193 // first: overall values for the whole SuperModule
194 Int_t iCalibMethod = 0;
195 Int_t iCalibPass = 0;
196 Float_t absoluteCalib = 0;
197 // third: info for each tower
198 Float_t relativeCalib[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
201 // just to make the initializations of the arrays are done correctly, let's use memset
202 memset(relativeCalib, 0, sizeof(relativeCalib));
204 // declare the branches
205 tree->SetBranchAddress("iSM", &iSM);
206 tree->SetBranchAddress("CalibMethod", &iCalibMethod);
207 tree->SetBranchAddress("CalibPass", &iCalibPass);
208 tree->SetBranchAddress("AbsoluteCalib", &absoluteCalib);
210 tree->SetBranchAddress("RelativeCalib", relativeCalib);
212 // indices for looping over the towers
216 for (int ient=0; ient<tree->GetEntries(); ient++) {
217 tree->GetEntry(ient);
219 // assume the index SuperModules come in order: i=iSM
220 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[iSM];
222 t->SetSuperModuleNum(iSM);
223 // first, overall values
224 t->SetCalibMethod(iCalibMethod);
225 t->SetCalibPass(iCalibPass);
226 t->SetAbsoluteCalib(absoluteCalib);
228 // third: info for each tower
229 for (Int_t j=0; j<nAPDPerSM; j++) {
230 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
231 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
233 // help variables: possibly modified or swapped indices
236 // assume that this info is already swapped and done for this basis?
238 // C side, oriented differently than A side: swap is requested
239 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
240 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
243 t->SetRelativeCalib(iColMod, iRowMod, relativeCalib[iCol][iRow]);
246 } // loop over entries
251 //____________________________________________________________________________
252 void AliEMCALCalibAbs::WriteRootCalibAbsInfo(const TString &rootFileName,
255 // write data to root file. ; coordinates given on SuperModule basis
256 TFile destFile(rootFileName, "recreate");
257 if (destFile.IsZombie()) {
262 TTree *tree = new TTree("tree","");
264 // variables for filling the TTree
265 Int_t iSM = 0; // SuperModule index
266 // list of values to be written
267 // first: overall values for the whole SuperModule
268 Int_t iCalibMethod = 0;
269 Int_t iCalibPass = 0;
270 Float_t absoluteCalib = 0;
271 // third: info for each tower
272 Float_t relativeCalib[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
275 // just to make the initializations of the arrays are done correctly, let's use memset
276 memset(relativeCalib, 0, sizeof(relativeCalib));
278 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
279 // for looping over towers
283 // declare the branches
285 tree->Branch("iSM", &iSM, "iSM/I");
286 tree->Branch("CalibMethod", &iCalibMethod, "CalibMethod/I");
287 tree->Branch("CalibPass", &iCalibPass, "CalibPass/I");
288 tree->Branch("AbsoluteCalib", &absoluteCalib, "AbsoluteCalib/F");
289 // third: info for each tower; see if a 2D array works OK or if we'll have to use 1D arrays instead
290 tree->Branch( "RelativeCalib", &relativeCalib, Form("RelativeCalib[%d][%d]/F", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
292 for (iSM = 0; iSM < fNSuperModule; iSM++) {
293 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[iSM];
295 iSM = t->GetSuperModuleNum();
296 // first, overall values
297 iCalibMethod = t->GetCalibMethod();
298 iCalibPass = t->GetCalibPass();
299 absoluteCalib = t->GetAbsoluteCalib();
301 // third: info for each tower
302 for (Int_t j=0; j<nAPDPerSM; j++) {
303 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
304 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
306 // help variables: possibly modified or swapped indices
309 // assume that this info is already swapped and done for this basis?
311 // C side, oriented differently than A side: swap is requested
312 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
313 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
316 relativeCalib[iColMod][iRowMod] = t->GetRelativeCalib(iCol, iRow);
328 //____________________________________________________________________________
329 AliEMCALCalibAbs::~AliEMCALCalibAbs()
331 fSuperModuleData.Delete();
334 //____________________________________________________________________________
335 AliEMCALSuperModuleCalibAbs * AliEMCALCalibAbs::GetSuperModuleCalibAbsNum(Int_t supModIndex)const
336 { // getter via index
337 for (int i=0; i<fNSuperModule; i++) {
338 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
339 if (t->GetSuperModuleNum() == supModIndex) {
344 // if we arrived here, then nothing was found.. just return a NULL pointer