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 outputFile << t->GetSuperModuleNum() << endl;
132 // info for each tower
133 for (Int_t j=0; j<nAPDPerSM; j++) {
134 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
135 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
137 tempCoeff = t->GetTC(iCol, iRow);
138 iSrc = t->GetSrc(iCol, iRow);
141 // C side, oriented differently than A side: swap is requested
142 iCol = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
143 iRow = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
146 outputFile << iCol << " " << iRow
148 << " " << iSrc << endl;
158 //____________________________________________________________________________
159 void AliEMCALCalibTempCoeff::ReadRootCalibTempCoeffInfo(const TString &rootFileName,
162 //Read data from root file. ; coordinates given on SuperModule basis
163 TFile inputFile(rootFileName, "read");
165 TTree *tree = (TTree*) inputFile.Get("tree");
167 ReadTreeCalibTempCoeffInfo(tree, swapSides);
174 //____________________________________________________________________________
175 void AliEMCALCalibTempCoeff::ReadTreeCalibTempCoeffInfo(TTree *tree,
178 // how many SuperModule's worth of info do we have?
179 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
180 fNSuperModule = tree->GetEntries();
182 Int_t iSM = 0; // SuperModule index
183 // list of values to be read
184 // info for each tower
185 Float_t tempCoeff[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
186 Int_t iSrc[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
189 // just to make the initializations of the arrays are done correctly, let's use memset
190 memset(tempCoeff, 0, sizeof(tempCoeff));
191 memset(iSrc, 0, sizeof(iSrc));
193 // declare the branches
194 tree->SetBranchAddress("iSM", &iSM);
195 tree->SetBranchAddress("TempCoeff", tempCoeff);
196 tree->SetBranchAddress("Src", iSrc);
198 // indices for looping over the towers
202 for (int ient=0; ient<tree->GetEntries(); ient++) {
203 tree->GetEntry(ient);
205 // assume the index SuperModules come in order: i=iSM
206 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[iSM];
208 t->SetSuperModuleNum(iSM);
210 // third: info for each tower
211 for (Int_t j=0; j<nAPDPerSM; j++) {
212 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
213 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
215 // help variables: possibly modified or swapped indices
218 // assume that this info is already swapped and done for this basis?
220 // C side, oriented differently than A side: swap is requested
221 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
222 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
225 t->SetTC(iColMod, iRowMod, tempCoeff[iCol][iRow]);
226 t->SetSrc(iColMod, iRowMod, iSrc[iCol][iRow]);
229 } // loop over entries
234 //____________________________________________________________________________
235 void AliEMCALCalibTempCoeff::WriteRootCalibTempCoeffInfo(const TString &rootFileName,
238 // write data to root file. ; coordinates given on SuperModule basis
239 TFile destFile(rootFileName, "recreate");
240 if (destFile.IsZombie()) {
245 TTree *tree = new TTree("tree","");
247 // variables for filling the TTree
248 Int_t iSM = 0; // SuperModule index
249 // list of values to be written
251 Float_t tempCoeff[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows];
252 Int_t iSrc[AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows]; // end - all values
254 // just to make the initializations of the arrays are done correctly, let's use memset
255 memset(tempCoeff, 0, sizeof(tempCoeff));
256 memset(iSrc, 0, sizeof(iSrc));
258 Int_t nAPDPerSM = AliEMCALGeoParams::fgkEMCALCols * AliEMCALGeoParams::fgkEMCALRows;
259 // for looping over towers
263 // declare the branches
265 tree->Branch("iSM", &iSM, "iSM/I");
266 // info for each tower; see if a 2D array works OK or if we'll have to use 1D arrays instead
267 tree->Branch( "TempCoeff", &tempCoeff, Form("TempCoeff[%d][%d]/F", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
268 tree->Branch( "Src", &iSrc, Form("Src[%d][%d]/I", AliEMCALGeoParams::fgkEMCALCols, AliEMCALGeoParams::fgkEMCALRows) );
270 for (iSM = 0; iSM < fNSuperModule; iSM++) {
271 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[iSM];
273 iSM = t->GetSuperModuleNum();
275 // info for each tower
276 for (Int_t j=0; j<nAPDPerSM; j++) {
277 iCol = j / AliEMCALGeoParams::fgkEMCALRows;
278 iRow = j % AliEMCALGeoParams::fgkEMCALRows;
280 // help variables: possibly modified or swapped indices
283 // assume that this info is already swapped and done for this basis?
285 // C side, oriented differently than A side: swap is requested
286 iColMod = AliEMCALGeoParams::fgkEMCALCols-1 - iCol;
287 iRowMod = AliEMCALGeoParams::fgkEMCALRows-1 - iRow;
290 tempCoeff[iColMod][iRowMod] = t->GetTC(iCol, iRow);
291 iSrc[iColMod][iRowMod] = t->GetSrc(iCol, iRow);
303 //____________________________________________________________________________
304 AliEMCALCalibTempCoeff::~AliEMCALCalibTempCoeff()
306 fSuperModuleData.Delete();
309 //____________________________________________________________________________
310 AliEMCALSuperModuleCalibTempCoeff * AliEMCALCalibTempCoeff::GetSuperModuleCalibTempCoeffNum(Int_t supModIndex)const
311 { // getter via index
312 for (int i=0; i<fNSuperModule; i++) {
313 AliEMCALSuperModuleCalibTempCoeff * t = (AliEMCALSuperModuleCalibTempCoeff*) fSuperModuleData[i];
314 if (t->GetSuperModuleNum() == supModIndex) {
319 // if we arrived here, then nothing was found.. just return a NULL pointer