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 temperature-dependence coefficients
26 #include "AliEMCALCalibTempCoeff.h"
30 ClassImp(AliEMCALCalibTempCoeff)
32 //____________________________________________________________________________
33 AliEMCALCalibTempCoeff::AliEMCALCalibTempCoeff(const int nSM) :
37 //Default constructor.
38 for (int i=0; i<fNSuperModule; i++) {
39 fSuperModuleData.Add(new AliEMCALSuperModuleCalibTempCoeff(i));
41 fSuperModuleData.Compress(); // compress the TObjArray
42 fSuperModuleData.SetOwner(kTRUE);
45 //____________________________________________________________________________
46 void AliEMCALCalibTempCoeff::ReadTextCalibTempCoeffInfo(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("AliEMCALCalibTempCoeff::ReadCalibTempCoeffInfo - Cannot open the APD info file %s\n", txtFileName.Data());
59 Int_t iSM = 0; // SuperModule index
63 // list of values to be read
65 Float_t tempCoeff = 0;
68 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
70 for (Int_t i = 0; i < fNSuperModule; i++) {
71 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[i];
73 printf("AliEMCALCalibTempCoeff::ReadCalibTempCoeffInfo - Error while reading input file; likely EOF..\n");
77 t->SetSuperModuleNum(iSM);
79 // info for each tower
80 for (Int_t j=0; j<nAPDPerSM; j++) {
81 inputFile >> iCol >> iRow >> tempCoeff >> iSrc;
83 // check that input values are not out bounds
84 if (iCol<0 || iCol>(AliEMCALGeoParams::fgkEMCALCols-1) ||
85 iRow<0 || iRow>(AliEMCALGeoParams::fgkEMCALRows-1) ) {
86 printf("AliEMCALCalibTempCoeff::ReadCalibTempCoeffInfo - Error while reading input file; j %d iCol %d iRow %d\n", j, iCol, iRow);
90 // assume that this info is already swapped and done for this basis?
92 // C side, oriented differently than A side: swap is requested
93 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
94 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
97 t->SetTC(iCol, iRow, tempCoeff);
98 t->SetSrc(iCol, iRow, iSrc);
108 //____________________________________________________________________________
109 void AliEMCALCalibTempCoeff::WriteTextCalibTempCoeffInfo(const TString &txtFileName,
112 // write data to txt file. ; coordinates given on SuperModule basis
114 std::ofstream outputFile(txtFileName.Data());
116 printf("AliEMCALCalibTempCoeff::WriteCalibTempCoeffInfo - Cannot open the APD output file %s\n", txtFileName.Data());
123 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
124 Float_t tempCoeff = 0;
127 for (Int_t i = 0; i < fNSuperModule; i++) {
128 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[i];
130 // info for each tower
131 for (Int_t j=0; j<nAPDPerSM; j++) {
132 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
133 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
135 tempCoeff = t->GetTC(iCol, iRow);
136 iSrc = t->GetSrc(iCol, iRow);
139 // C side, oriented differently than A side: swap is requested
140 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
141 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
144 outputFile << iCol << " " << iRow
146 << " " << iSrc << endl;
156 //____________________________________________________________________________
157 void AliEMCALCalibTempCoeff::ReadRootCalibTempCoeffInfo(const TString &rootFileName,
160 //Read data from root file. ; coordinates given on SuperModule basis
161 TFile inputFile(rootFileName, "read");
163 TTree *tree = (TTree*) inputFile.Get("tree");
165 ReadTreeCalibTempCoeffInfo(tree, swapSides);
172 //____________________________________________________________________________
173 void AliEMCALCalibTempCoeff::ReadTreeCalibTempCoeffInfo(TTree *tree,
176 // how many SuperModule's worth of info do we have?
177 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
178 fNSuperModule = tree->GetEntries();
180 Int_t iSM = 0; // SuperModule index
181 // list of values to be read
182 // info for each tower
183 Float_t tempCoeff[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
184 Int_t iSrc[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
187 // just to make the initializations of the arrays are done correctly, let's use memset
188 memset(tempCoeff, 0, sizeof(tempCoeff));
189 memset(iSrc, 0, sizeof(iSrc));
191 // declare the branches
192 tree->SetBranchAddress("iSM", &iSM);
193 tree->SetBranchAddress("TempCoeff", tempCoeff);
194 tree->SetBranchAddress("Src", iSrc);
196 // indices for looping over the towers
200 for (int ient=0; ient<tree->GetEntries(); ient++) {
201 tree->GetEntry(ient);
203 // assume the index SuperModules come in order: i=iSM
204 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[iSM];
206 t->SetSuperModuleNum(iSM);
208 // third: info for each tower
209 for (Int_t j=0; j<nAPDPerSM; j++) {
210 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
211 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
213 // help variables: possibly modified or swapped indices
216 // assume that this info is already swapped and done for this basis?
218 // C side, oriented differently than A side: swap is requested
219 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
220 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
223 t->SetTC(iColMod, iRowMod, tempCoeff[iCol][iRow]);
224 t->SetSrc(iColMod, iRowMod, iSrc[iCol][iRow]);
227 } // loop over entries
232 //____________________________________________________________________________
233 void AliEMCALCalibTempCoeff::WriteRootCalibTempCoeffInfo(const TString &rootFileName,
236 // write data to root file. ; coordinates given on SuperModule basis
237 TFile destFile(rootFileName, "recreate");
238 if (destFile.IsZombie()) {
243 TTree *tree = new TTree("tree","");
245 // variables for filling the TTree
246 Int_t iSM = 0; // SuperModule index
247 // list of values to be written
249 Float_t tempCoeff[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
250 Int_t iSrc[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows]; // end - all values
252 // just to make the initializations of the arrays are done correctly, let's use memset
253 memset(tempCoeff, 0, sizeof(tempCoeff));
254 memset(iSrc, 0, sizeof(iSrc));
256 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
257 // for looping over towers
261 // declare the branches
263 tree->Branch("iSM", &iSM, "iSM/I");
264 // info for each tower; see if a 2D array works OK or if we'll have to use 1D arrays instead
265 tree->Branch( "TempCoeff", &tempCoeff, Form("TempCoeff[%d][%d]/F", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
266 tree->Branch( "Src", &iSrc, Form("Src[%d][%d]/I", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
268 for (iSM = 0; iSM < fNSuperModule; iSM++) {
269 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[iSM];
271 iSM = t->GetSuperModuleNum();
273 // info for each tower
274 for (Int_t j=0; j<nAPDPerSM; j++) {
275 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
276 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
278 // help variables: possibly modified or swapped indices
281 // assume that this info is already swapped and done for this basis?
283 // C side, oriented differently than A side: swap is requested
284 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
285 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
288 tempCoeff[iColMod][iRowMod] = t->GetTC(iCol, iRow);
289 iSrc[iColMod][iRowMod] = t->GetSrc(iCol, iRow);
301 //____________________________________________________________________________
302 AliEMCALCalibTempCoeff::~AliEMCALCalibTempCoeff()
304 fSuperModuleData.Delete();
307 //____________________________________________________________________________
308 AliEMCALSuperModuleCalibTempCoeff * AliEMCALCalibTempCoeff::GetSuperModuleCalibTempCoeffNum(Int_t supModIndex)const
310 for (int i=0; i<fNSuperModule; i++) {
311 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[i];
312 if (t->GetSuperModuleNum() == supModIndex) {
317 // if we arrived here, then nothing was found.. just return a NULL pointer