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
44 //____________________________________________________________________________
45 void AliEMCALCalibAbs::ReadTextCalibAbsInfo(Int_t nSM, const TString &txtFileName,
48 //Read data from txt file. ; coordinates given on SuperModule basis
50 std::ifstream inputFile(txtFileName.Data());
52 printf("AliEMCALCalibAbs::ReadCalibAbsInfo - Cannot open the APD info file %s\n", txtFileName.Data());
58 Int_t iSM = 0; // SuperModule index
62 // list of values to be read
63 // first: overall values for the whole SuperModule
66 Float_t absoluteCalib;
67 // third: info for each tower
68 Float_t relativeCalib; // (ADC>GeV relative gain/conversion), value around 1
71 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
73 for (Int_t i = 0; i < fNSuperModule; i++) {
74 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
76 printf("AliEMCALCalibAbs::ReadCalibAbsInfo - Error while reading input file; likely EOF..");
80 t->SetSuperModuleNum(iSM);
82 // first: overall values for the whole SuperModule
83 inputFile >> iCalibMethod >> iCalibPass >> absoluteCalib;
84 t->SetCalibMethod(iCalibMethod);
85 t->SetCalibPass(iCalibPass);
86 t->SetAbsoluteCalib(absoluteCalib);
88 // third: info for each tower
89 for (Int_t j=0; j<nAPDPerSM; j++) {
90 inputFile >> iCol >> iRow
93 // assume that this info is already swapped and done for this basis?
95 // C side, oriented differently than A side: swap is requested
96 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
97 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
100 t->SetRelativeCalib(iCol, iRow, relativeCalib);
110 //____________________________________________________________________________
111 void AliEMCALCalibAbs::WriteTextCalibAbsInfo(const TString &txtFileName,
114 // write data to txt file. ; coordinates given on SuperModule basis
116 std::ofstream outputFile(txtFileName.Data());
118 printf("AliEMCALCalibAbs::WriteCalibAbsInfo - Cannot open the APD output file %s\n", txtFileName.Data());
125 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
126 Float_t relativeCalib = 0;
127 for (Int_t i = 0; i < fNSuperModule; i++) {
128 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
130 // first: overall values for the whole SuperModule
131 outputFile << t->GetSuperModuleNum() << endl;
132 outputFile << t->GetCalibMethod() << " "
133 << t->GetCalibPass() << " "
134 << t->GetAbsoluteCalib() << endl;
136 // third: info for each tower
137 for (Int_t j=0; j<nAPDPerSM; j++) {
138 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
139 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
141 relativeCalib = t->GetRelativeCalib(iCol, iRow);
144 // C side, oriented differently than A side: swap is requested
145 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
146 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
149 outputFile << iCol << " " << iRow
150 << " " << relativeCalib << endl;
160 //____________________________________________________________________________
161 void AliEMCALCalibAbs::ReadRootCalibAbsInfo(const TString &rootFileName,
164 //Read data from root file. ; coordinates given on SuperModule basis
165 TFile inputFile(rootFileName, "read");
167 TTree *tree = (TTree*) inputFile.Get("tree");
169 ReadTreeCalibAbsInfo(tree, swapSides);
176 //____________________________________________________________________________
177 void AliEMCALCalibAbs::ReadTreeCalibAbsInfo(TTree *tree,
180 // how many SuperModule's worth of info do we have?
181 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
182 fNSuperModule = tree->GetEntries();
184 Int_t iSM = 0; // SuperModule index
185 // list of values to be read
186 // first: overall values for the whole SuperModule
188 Int_t iCalibPass = 0;
189 Float_t absoluteCalib = 0;
190 // third: info for each tower
191 Float_t relativeCalib[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
194 // just to make the initializations of the arrays are done correctly, let's use memset
195 memset(relativeCalib, 0, sizeof(relativeCalib));
197 // declare the branches
198 tree->SetBranchAddress("iSM", &iSM);
199 tree->SetBranchAddress("CalibMethod", &iCalibMethod);
200 tree->SetBranchAddress("CalibPass", &iCalibPass);
201 tree->SetBranchAddress("AbsoluteCalib", &absoluteCalib);
203 tree->SetBranchAddress("RelativeCalib", relativeCalib);
205 // indices for looping over the towers
209 for (int ient=0; ient<tree->GetEntries(); ient++) {
210 tree->GetEntry(ient);
212 // assume the index SuperModules come in order: i=iSM
213 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[iSM];
215 t->SetSuperModuleNum(iSM);
216 // first, overall values
217 t->SetCalibMethod(iCalibMethod);
218 t->SetCalibPass(iCalibPass);
219 t->SetAbsoluteCalib(absoluteCalib);
221 // third: info for each tower
222 for (Int_t j=0; j<nAPDPerSM; j++) {
223 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
224 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
226 // help variables: possibly modified or swapped indices
229 // assume that this info is already swapped and done for this basis?
231 // C side, oriented differently than A side: swap is requested
232 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
233 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
236 t->SetRelativeCalib(iColMod, iRowMod, relativeCalib[iCol][iRow]);
239 } // loop over entries
244 //____________________________________________________________________________
245 void AliEMCALCalibAbs::WriteRootCalibAbsInfo(const TString &rootFileName,
248 // write data to root file. ; coordinates given on SuperModule basis
249 TFile destFile(rootFileName, "recreate");
250 if (destFile.IsZombie()) {
255 TTree *tree = new TTree("tree","");
257 // variables for filling the TTree
258 Int_t iSM = 0; // SuperModule index
259 // list of values to be written
260 // first: overall values for the whole SuperModule
261 Int_t iCalibMethod = 0;
262 Int_t iCalibPass = 0;
263 Float_t absoluteCalib = 0;
264 // third: info for each tower
265 Float_t relativeCalib[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
268 // just to make the initializations of the arrays are done correctly, let's use memset
269 memset(relativeCalib, 0, sizeof(relativeCalib));
271 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
272 // for looping over towers
276 // declare the branches
278 tree->Branch("iSM", &iSM, "iSM/I");
279 tree->Branch("CalibMethod", &iCalibMethod, "CalibMethod/I");
280 tree->Branch("CalibPass", &iCalibPass, "CalibPass/I");
281 tree->Branch("AbsoluteCalib", &absoluteCalib, "AbsoluteCalib/F");
282 // third: info for each tower; see if a 2D array works OK or if we'll have to use 1D arrays instead
283 tree->Branch( "RelativeCalib", &relativeCalib, Form("RelativeCalib[%d][%d]/F", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
285 for (iSM = 0; iSM < fNSuperModule; iSM++) {
286 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[iSM];
288 iSM = t->GetSuperModuleNum();
289 // first, overall values
290 iCalibMethod = t->GetCalibMethod();
291 iCalibPass = t->GetCalibPass();
292 absoluteCalib = t->GetAbsoluteCalib();
294 // third: info for each tower
295 for (Int_t j=0; j<nAPDPerSM; j++) {
296 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
297 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
299 // help variables: possibly modified or swapped indices
302 // assume that this info is already swapped and done for this basis?
304 // C side, oriented differently than A side: swap is requested
305 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
306 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
309 relativeCalib[iColMod][iRowMod] = t->GetRelativeCalib(iCol, iRow);
321 //____________________________________________________________________________
322 AliEMCALCalibAbs::~AliEMCALCalibAbs()
324 fSuperModuleData.Delete();
327 //____________________________________________________________________________
328 AliEMCALSuperModuleCalibAbs * AliEMCALCalibAbs::GetSuperModuleCalibAbsNum(Int_t supModIndex)const
330 for (int i=0; i<fNSuperModule; i++) {
331 AliEMCALSuperModuleCalibAbs * t = (AliEMCALSuperModuleCalibAbs*) fSuperModuleData[i];
332 if (t->GetSuperModuleNum() == supModIndex) {
337 // if we arrived here, then nothing was found.. just return a NULL pointer