EcalTBTools/src/TB06TreeH2.cc

0001 #include "Calibration/EcalTBTools/interface/TB06RecoH2.h"
0002 #include "Calibration/EcalTBTools/interface/TB06TreeH2.h"
0003 #include "TFile.h"
0004 #include "TTree.h"
0005 #include "TClonesArray.h"
0006
0007 #include <iostream>
0008
0009 TB06TreeH2::TB06TreeH2(const std::string &fileName, const std::string &treeName)
0010     : m_file(nullptr), m_tree(nullptr), m_data(nullptr), m_dataSize(0) {
0011   TDirectory *dir = gDirectory;
0012   m_file = new TFile(fileName.c_str(), "RECREATE");
0013   m_file->cd();
0014   m_tree = new TTree(treeName.c_str(), "Analysis tree");
0015   m_tree->SetAutoSave(10000000);
0016   dir->cd();
0017
0018   //  m_tree->cd () ;
0019   m_data = new TClonesArray(TB06RecoH2::Class(), 1);
0020   m_data->ExpandCreateFast(1);
0021
0022   //  m_tree->Branch ("EGCO", &m_data, 64000, 2) ;
0023   m_tree->Branch("TB06O", &m_data, 64000, 2);
0024   m_tree->Print();
0025 }
0026
0027 // -------------------------------------------------------------------
0028
0029 TB06TreeH2::~TB06TreeH2() {
0030   std::cout << "[TB06TreeH2][dtor] saving TTree " << m_tree->GetName() << " with " << m_tree->GetEntries() << " entries"
0031             << " on file: " << m_file->GetName() << std::endl;
0032
0033   m_file->Write();
0034   delete m_tree;
0035   m_file->Close();
0036   delete m_file;
0037   delete m_data;
0038 }
0039
0040 // -------------------------------------------------------------------
0041
0042 //! to be called at each loop
0043 void TB06TreeH2::store(const int &tableIsMoving,
0044                        const int &run,
0045                        const int &event,
0046                        const int &S6adc,
0047                        const double &xhodo,
0048                        const double &yhodo,
0049                        const double &xslope,
0050                        const double &yslope,
0051                        const double &xquality,
0052                        const double &yquality,
0053                        const int &icMax,
0054                        const int &ietaMax,
0055                        const int &iphiMax,
0056                        const double &beamEnergy,
0057                        const double ampl[49],
0058                        const int &wcAXo,
0059                        const int &wcAYo,
0060                        const int &wcBXo,
0061                        const int &wcBYo,
0062                        const int &wcCXo,
0063                        const int &wcCYo,
0064                        const double &xwA,
0065                        const double &ywA,
0066                        const double &xwB,
0067                        const double &ywB,
0068                        const double &xwC,
0069                        const double &ywC,
0070                        const float &S1adc,
0071                        const float &S2adc,
0072                        const float &S3adc,
0073                        const float &S4adc,
0074                        const float &VM1,
0075                        const float &VM2,
0076                        const float &VM3,
0077                        const float &VM4,
0078                        const float &VM5,
0079                        const float &VM6,
0080                        const float &VM7,
0081                        const float &VM8,
0082                        const float &VMF,
0083                        const float &VMB,
0084                        const float &CK1,
0085                        const float &CK2,
0086                        const float &CK3,
0087                        const float &BH1,
0088                        const float &BH2,
0089                        const float &BH3,
0090                        const float &BH4,
0091                        const float &TOF1S,
0092                        const float &TOF2S,
0093                        const float &TOF1J,
0094                        const float &TOF2J) {
0095   m_data->Clear();
0096   TB06RecoH2 *entry = static_cast<TB06RecoH2 *>(m_data->AddrAt(0));
0097
0098   entry->reset();
0099   //  reset (entry->myCalibrationMap) ;
0100
0101   entry->tableIsMoving = tableIsMoving;
0102   entry->run = run;
0103   entry->event = event;
0104   entry->S6ADC = S6adc;
0105
0106   entry->MEXTLindex = icMax;
0107   entry->MEXTLeta = ietaMax;
0108   entry->MEXTLphi = iphiMax;
0109   entry->MEXTLenergy = ampl[24];
0110   entry->beamEnergy = beamEnergy;
0111
0112   for (int eta = 0; eta < 7; ++eta)
0113     for (int phi = 0; phi < 7; ++phi) {
0114       // FIXME capire l'orientamento di phi!
0115       // FIXME capire se eta, phi iniziano da 1 o da 0
0116       entry->localMap[eta][phi] = ampl[eta * 7 + phi];
0117     }
0118
0119   //[Edgar] S1 uncleaned, uncalibrated energy
0120   entry->S1uncalib_ = ampl[24];
0121
0122   //[Edgar] S25 uncleaned, uncalibrated energy
0123   for (int eta = 1; eta < 6; ++eta)
0124     for (int phi = 1; phi < 6; ++phi) {
0125       entry->S25uncalib_ += entry->localMap[eta][phi];
0126     }
0127
0128   //[Edgar] S49 uncleaned, uncalibrated energy
0129   for (int eta = 0; eta < 7; ++eta)
0130     for (int phi = 0; phi < 7; ++phi) {
0131       entry->S49uncalib_ += entry->localMap[eta][phi];
0132     }
0133
0134   //[Edgar] S9 uncleaned, uncalibrated energy
0135   for (int eta = 2; eta < 5; ++eta)
0136     for (int phi = 2; phi < 5; ++phi) {
0137       entry->S9uncalib_ += entry->localMap[eta][phi];
0138     }
0139
0140   entry->xHodo = xhodo;
0141   entry->yHodo = yhodo;
0142   entry->xSlopeHodo = xslope;
0143   entry->ySlopeHodo = yslope;
0144   entry->xQualityHodo = xquality;
0145   entry->yQualityHodo = yquality;
0146   entry->wcAXo_ = wcAXo;
0147   entry->wcAYo_ = wcAYo;
0148   entry->wcBXo_ = wcBXo;
0149   entry->wcBYo_ = wcBYo;
0150   entry->wcCXo_ = wcCXo;
0151   entry->wcCYo_ = wcCYo;
0152   entry->xwA_ = xwA;
0153   entry->ywA_ = ywA;
0154   entry->xwB_ = xwB;
0155   entry->ywB_ = ywB;
0156   entry->xwC_ = xwC;
0157   entry->ywC_ = ywC;
0158   entry->S1adc_ = S1adc;
0159   entry->S2adc_ = S2adc;
0160   entry->S3adc_ = S3adc;
0161   entry->S4adc_ = S4adc;
0162   entry->VM1_ = VM1;
0163   entry->VM2_ = VM2;
0164   entry->VM3_ = VM3;
0165   entry->VM4_ = VM4;
0166   entry->VM5_ = VM5;
0167   entry->VM6_ = VM6;
0168   entry->VM7_ = VM7;
0169   entry->VM8_ = VM8;
0170   entry->VMF_ = VMF;
0171   entry->VMB_ = VMB;
0172   entry->CK1_ = CK1;
0173   entry->CK2_ = CK2;
0174   entry->CK3_ = CK3;
0175   entry->BH1_ = BH1;
0176   entry->BH2_ = BH2;
0177   entry->BH3_ = BH3;
0178   entry->BH4_ = BH4;
0179   entry->TOF1S_ = TOF1S;
0180   entry->TOF2S_ = TOF2S;
0181   entry->TOF1J_ = TOF1J;
0182   entry->TOF2J_ = TOF2J;
0183
0184   entry->convFactor = 0.;
0185
0186   /*
0187   // loop over the 5x5 see (1)
0188   for (int xtal=0 ; xtal<25 ; ++xtal)
0189     {
0190       int ieta = xtal/5 + 3 ;
0191       int iphi = xtal%5 + 8 ;
0192       entry->myCalibrationMap[ieta][iphi] = ampl[xtal] ;
0193     } // loop over the 5x5
0194
0195   entry->electron_Tr_Pmag_ = beamEnergy ;
0196
0197         entry->centralCrystalEta_ = ietaMax ;
0198         entry->centralCrystalPhi_ = iphiMax ;
0199         entry->centralCrystalEnergy_ = ampl[12] ;
0200
0201   // this is a trick
0202   entry->electron_Tr_Peta_ = xhodo ;
0203   entry->electron_Tr_Pphi_ = yhodo ;
0204   */
0205   m_tree->Fill();
0206 }
0207
0208 // -------------------------------------------------------------------
0209
0210 void TB06TreeH2::reset(float crystal[11][21]) {
0211   for (int eta = 0; eta < 11; ++eta) {
0212     for (int phi = 0; phi < 21; ++phi) {
0213       crystal[eta][phi] = -999.;
0214     }
0215   }
0216 }
0217
0218 // -------------------------------------------------------------------
0219
0220 void TB06TreeH2::check() {
0221   TB06RecoH2 *entry = static_cast<TB06RecoH2 *>(m_data->AddrAt(0));
0222
0223   std::cout << "[TB06TreeH2][check]reading . . . \n";
0224   std::cout << "[TB06TreeH2][check] entry->run: " << entry->run << "\n";
0225   std::cout << "[TB06TreeH2][check] entry->event: " << entry->event << "\n";
0226   std::cout << "[TB06TreeH2][check] entry->tableIsMoving: " << entry->tableIsMoving << "\n";
0227   std::cout << "[TB06TreeH2][check] entry->MEXTLeta: " << entry->MEXTLeta << "\n";
0228   std::cout << "[TB06TreeH2][check] entry->MEXTLphi: " << entry->MEXTLphi << "\n";
0229   std::cout << "[TB06TreeH2][check] entry->MEXTLenergy: " << entry->MEXTLenergy << "\n";
0230
0231   for (int eta = 0; eta < 7; ++eta)
0232     for (int phi = 0; phi < 7; ++phi)
0233       std::cout << "[TB06TreeH2][check]   entry->localMap[" << eta << "][" << phi << "]: " << entry->localMap[eta][phi]
0234                 << "\n";
0235
0236   std::cout << "[TB06TreeH2][check] entry->xHodo: " << entry->xHodo << "\n";
0237   std::cout << "[TB06TreeH2][check] entry->yHodo: " << entry->yHodo << "\n";
0238   std::cout << "[TB06TreeH2][check] entry->xSlopeHodo: " << entry->xSlopeHodo << "\n";
0239   std::cout << "[TB06TreeH2][check] entry->ySlopeHodo: " << entry->ySlopeHodo << "\n";
0240   std::cout << "[TB06TreeH2][check] entry->xQualityHodo: " << entry->xQualityHodo << "\n";
0241   std::cout << "[TB06TreeH2][check] entry->yQualityHodo: " << entry->yQualityHodo << "\n";
0242   std::cout << "[TB06TreeH2][check] entry->convFactor: " << entry->convFactor << "\n";
0243
0244   /* to be implemented with the right variables
0245   std::cout << "[TB06TreeH2][check] ------------------------" << std::endl ;
0246   std::cout << "[TB06TreeH2][check] " << entry->variable_name << std::endl ;
0247   */
0248 }
0249
0250 /* (1) to fill the 25 crystals vector
0251
0252    for (UInt_t icry=0 ; icry<25 ; ++icry)
0253      {
0254        UInt_t row = icry / 5 ;
0255        Int_t column = icry % 5 ;
0256        try
0257            {
0258              EBDetId tempo (maxHitId.ieta()+column-2,
0259                             maxHitId.iphi()+row-2,
0260                             EBDetId::ETAPHIMODE) ;
0261
0262              Xtals5x5.push_back (tempo) ;
0263              amplitude [icry] = hits->find (Xtals5x5[icry])->energy () ;
0264
0265            }
0266        catch ( std::runtime_error &e )
0267            {
0268              std::cout << "Cannot construct 5x5 matrix around EBDetId "
0269                      << maxHitId << std::endl ;
0270              return ;
0271            }
0272      } // loop over the 5x5 matrix
0273
0274
0275 */