Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:23:13

0001 /*******************************************
0002 
0003  DumpLaserDB.cpp
0004  Author: Giovanni.Organtini@roma1.infn.it
0005  Date:   may 2011
0006 
0007  Description: reads the laser DB and builds
0008    a plain root ntuple with primitives and
0009    associated quantities.
0010 
0011  ******************************************/
0012 #include <climits>
0013 #include <cstdlib>
0014 #include <ctime>
0015 #include <iomanip>
0016 #include <iostream>
0017 #include <sstream>
0018 #include <string>
0019 #include <vector>
0020 
0021 #include <TFile.h>
0022 #include <TTree.h>
0023 
0024 #include "boost/program_options.hpp"
0025 
0026 #include "OnlineDB/EcalCondDB/interface/EcalCondDBInterface.h"
0027 #include "OnlineDB/EcalCondDB/interface/all_lmf_types.h"
0028 
0029 typedef struct lmffdata {
0030   int detId;
0031   int logic_id;
0032   int eta;
0033   int phi;
0034   int x;
0035   int y;
0036   int z;
0037   int nevts;
0038   int quality;
0039   float pnMean;
0040   float pnRMS;
0041   float pnM3;
0042   float pnAoverpnBMean;
0043   float pnAoverpnBRMS;
0044   float pnAoverpnBM3;
0045   int pnFlag;
0046   float Mean;
0047   float RMS;
0048   float M3;
0049   float APDoverPNAMean;
0050   float APDoverPNARMS;
0051   float APDoverPNAM3;
0052   float APDoverPNBMean;
0053   float APDoverPNBRMS;
0054   float APDoverPNBM3;
0055   float APDoverPNMean;
0056   float APDoverPNRMS;
0057   float APDoverPNM3;
0058   float alpha;
0059   float beta;
0060   int flag;
0061 } LMFData;
0062 
0063 typedef struct matacq {
0064   int fit_method;
0065   float mtq_ampl;
0066   float mtq_time;
0067   float mtq_rise;
0068   float mtq_fwhm;
0069   float mtq_fw20;
0070   float mtq_fw80;
0071   float mtq_sliding;
0072 } LMFMtq;
0073 
0074 typedef struct laser_config {
0075   int wavelength;
0076   float vfe_gain;
0077   float pn_gain;
0078   float lsr_power;
0079   float lsr_attenuator;
0080   float lsr_current;
0081   float lsr_delay_1;
0082   float lsr_delay_2;
0083 } LMFLaserConfig;
0084 
0085 typedef struct DetIdData {
0086   int hashedIndex;
0087   int eta;
0088   int phi;
0089   int x;
0090   int y;
0091   int z;
0092 } detIdData;
0093 
0094 class CondDBApp {
0095 public:
0096   /**
0097    *   App constructor; Makes the database connection
0098    */
0099   CondDBApp(std::string sid, std::string user, std::string pass, run_t r1, run_t r2) {
0100     try {
0101       std::cout << "Making connection to " << sid << " using username " << user << std::flush;
0102       econn = new EcalCondDBInterface(sid, user, pass);
0103       std::cout << "Done." << std::endl;
0104       std::cout << "Getting data for run" << std::flush;
0105       if (r2 <= 0) {
0106         std::cout << " " << r1 << std::endl;
0107       } else {
0108         std::cout << "s between " << r1 << " and " << r2 << std::endl;
0109       }
0110       run_min = r1;
0111       run_max = r2;
0112     } catch (std::runtime_error &e) {
0113       std::cerr << e.what() << std::endl;
0114       exit(-1);
0115     }
0116   }
0117 
0118   /**
0119    *  App destructor;  Cleans up database connection
0120    */
0121   ~CondDBApp() {
0122     tree->Write();
0123     tfile->Close();
0124     delete econn;
0125   }
0126 
0127   void zero(std::vector<LMFMtq> &lmfmtq) {
0128     for (unsigned int i = 0; i < lmfmtq.size(); i++) {
0129       lmfmtq[i].fit_method = -1;
0130       lmfmtq[i].mtq_ampl = -1000;
0131       lmfmtq[i].mtq_time = -1000;
0132       lmfmtq[i].mtq_rise = -1000;
0133       lmfmtq[i].mtq_fwhm = -1000;
0134       lmfmtq[i].mtq_fw20 = -1000;
0135       lmfmtq[i].mtq_fw80 = -1000;
0136       lmfmtq[i].mtq_sliding = -1000;
0137     }
0138   }
0139 
0140   void zero(std::vector<LMFLaserConfig> &lmfConf) {
0141     for (unsigned int i = 0; i < lmfConf.size(); i++) {
0142       lmfConf[i].wavelength = -1;
0143       lmfConf[i].vfe_gain = -1000;
0144       lmfConf[i].pn_gain = -1000;
0145       lmfConf[i].lsr_power = -1000;
0146       lmfConf[i].lsr_attenuator = -1000;
0147       lmfConf[i].lsr_current = -1000;
0148       lmfConf[i].lsr_delay_1 = -1000;
0149       lmfConf[i].lsr_delay_2 = -1000;
0150     }
0151   }
0152 
0153   void zero(std::vector<LMFData> &lmfdata) {
0154     for (unsigned int i = 0; i < lmfdata.size(); i++) {
0155       lmfdata[i].detId = -1;
0156       lmfdata[i].logic_id = -1;
0157       lmfdata[i].eta = -1000;
0158       lmfdata[i].phi = -1000;
0159       lmfdata[i].x = -1000;
0160       lmfdata[i].y = -1000;
0161       lmfdata[i].z = -1000;
0162       lmfdata[i].nevts = 0;
0163       lmfdata[i].quality = -1;
0164       lmfdata[i].pnMean = -1000.;
0165       lmfdata[i].pnRMS = -1000.;
0166       lmfdata[i].pnM3 = -1000.;
0167       lmfdata[i].pnAoverpnBMean = -1000.;
0168       lmfdata[i].pnAoverpnBRMS = -1000.;
0169       lmfdata[i].pnAoverpnBM3 = -1000.;
0170       lmfdata[i].pnFlag = -1000;
0171       lmfdata[i].Mean = -1000.;
0172       lmfdata[i].RMS = -1000.;
0173       lmfdata[i].M3 = -1000.;
0174       lmfdata[i].APDoverPNAMean = -1000.;
0175       lmfdata[i].APDoverPNARMS = -1000.;
0176       lmfdata[i].APDoverPNAM3 = -1000.;
0177       lmfdata[i].APDoverPNBMean = -1000.;
0178       lmfdata[i].APDoverPNBRMS = -1000.;
0179       lmfdata[i].APDoverPNBM3 = -1000.;
0180       lmfdata[i].APDoverPNMean = -1000.;
0181       lmfdata[i].APDoverPNRMS = -1000.;
0182       lmfdata[i].APDoverPNM3 = -1000.;
0183       lmfdata[i].alpha = -1000.;
0184       lmfdata[i].beta = -1000.;
0185       lmfdata[i].flag = -1000;
0186     }
0187   }
0188 
0189   detIdData getCoords(int detid) {
0190     detIdData ret;
0191     DetId d(detid);
0192     if (d.subdetId() == EcalBarrel) {
0193       EBDetId eb(detid);
0194       ret.hashedIndex = eb.hashedIndex();
0195       ret.eta = eb.ieta();
0196       ret.phi = eb.iphi();
0197       ret.x = -1000;
0198       ret.y = -1000;
0199       ret.z = -1000;
0200     } else {
0201       EEDetId ee(detid);
0202       ret.hashedIndex = ee.hashedIndex() + 61200;
0203       ret.x = ee.ix();
0204       ret.y = ee.iy();
0205       ret.z = ee.zside();
0206       ret.eta = -1000;
0207       ret.phi = -1000;
0208     }
0209     return ret;
0210   }
0211 
0212   bool userCheck(LMFData d) {
0213     /*** EXPERTS ONLY ***/
0214     bool ret = true;
0215     return ret;
0216   }
0217 
0218   void doRun() {
0219     init();
0220     std::map<int, int> detids = econn->getLogicId2DetIdMap();
0221     std::cout << "Crystal map got: " << detids.size() << std::endl;
0222     for (run_t run = run_min; run <= run_max; run++) {
0223       std::cout << "Analyzing run " << run << std::endl << std::flush;
0224       run_num = run;
0225       // for each run collect
0226       LMFSeqDat s(econn);
0227       std::map<int, LMFSeqDat> sequences = s.fetchByRunNumber(run);
0228       std::map<int, LMFSeqDat>::const_iterator si = sequences.begin();
0229       std::map<int, LMFSeqDat>::const_iterator se = sequences.end();
0230       while (si != se) {
0231         // seq start
0232         seqStart = si->second.getSequenceStart().epoch();
0233         // seq stop
0234         seqStop = si->second.getSequenceStop().epoch();
0235         // seq number
0236         seqNum = si->second.getSequenceNumber();
0237         std::cout << std::endl << " Seq. " << seqNum;
0238         LMFRunIOV riov(econn);
0239         std::list<LMFRunIOV> run_iovs = riov.fetchBySequence(si->second);
0240         std::list<LMFRunIOV>::const_iterator ri = run_iovs.begin();
0241         std::list<LMFRunIOV>::const_iterator re = run_iovs.end();
0242         while (ri != re) {
0243           // lmr
0244           lmr = ri->getLmr();
0245           std::cout << std::endl << "  LMR " << std::setw(2) << lmr << std::flush;
0246           // color
0247           sprintf(&color[0], "%s", ri->getColorShortName().c_str());
0248           // subrun times
0249           subrunstart = ri->getSubRunStart().epoch();
0250           subrunstop = ri->getSubRunEnd().epoch();
0251           // blue laser configuration for this run
0252           LMFLaserPulseDat mtqConf(econn, "BLUE");
0253           mtqConf.setLMFRunIOV(*ri);
0254           mtqConf.fetch();
0255           std::list<int> mtqLogicIds = mtqConf.getLogicIds();
0256           LMFLaserConfigDat laserConfig(econn);
0257           laserConfig.setLMFRunIOV(*ri);
0258           laserConfig.fetch();
0259           std::list<int> laserConfigLogicIds = laserConfig.getLogicIds();
0260           // *** get data ***
0261           LMFPrimDat prim(econn, color, "LASER");
0262           prim.setLMFRunIOV(*ri);
0263           std::vector<std::string> channels;
0264           // uncomment the following line to selects just three channels
0265           /*
0266       channels.push_back("2012043034");
0267       channels.push_back("2012060033");
0268       channels.push_back("2012034040");
0269       prim.setWhereClause("(LOGIC_ID = :I1 OR LOGIC_ID = :I2 OR LOGIC_ID = :I3)", channels); // selects only endcap primitives
0270       */
0271           prim.fetch();
0272           if (!prim.getLogicIds().empty()) {
0273             LMFRunDat run_dat(econn);
0274             run_dat.setLMFRunIOV(*ri);
0275             /* uncomment the following to select only endcaps
0276         run_dat.setWhereClause("LOGIC_ID > 2000000000"); // selects only endcap primitives
0277         */
0278             run_dat.fetch();
0279             LMFPnPrimDat pnPrim(econn, color, "LASER");
0280             pnPrim.setLMFRunIOV(*ri);
0281             /* uncomment the following to select only endcaps
0282         pnPrim.setWhereClause("LOGIC_ID > 2000000000"); // selects only endcap primitives
0283         */
0284             pnPrim.fetch();
0285             // *** run dat ***
0286             std::list<int> logic_ids = run_dat.getLogicIds();
0287             std::list<int>::const_iterator li = logic_ids.begin();
0288             std::list<int>::const_iterator le = logic_ids.end();
0289             int count = 0;
0290             int xcount = 0;
0291             std::vector<LMFData> lmfdata(detids.size());
0292             std::vector<LMFMtq> lmfmtq(detids.size());
0293             std::vector<LMFLaserConfig> lmflasConf(detids.size());
0294             zero(lmfdata);
0295             zero(lmfmtq);
0296             zero(lmflasConf);
0297             while (li != le) {
0298               // here the logic_id's are those of a LMR: transform into crystals
0299               int logic_id = *li;
0300               std::vector<EcalLogicID> xtals = econn->getEcalLogicIDForLMR(logic_id);
0301               std::cout << " Size = " << std::setw(4) << xtals.size();
0302               for (unsigned int j = 0; j < xtals.size(); j++) {
0303                 int xtal_id = xtals[j].getLogicID();
0304                 int detId = detids[xtal_id];
0305                 detIdData dd = getCoords(detId);
0306                 int index = dd.hashedIndex;
0307                 lmfdata[index].nevts = run_dat.getEvents(logic_id);
0308                 lmfdata[index].detId = detId;
0309                 lmfdata[index].logic_id = xtal_id;
0310                 lmfdata[index].quality = run_dat.getQualityFlag(logic_id);
0311                 lmfdata[index].eta = dd.eta;
0312                 lmfdata[index].phi = dd.phi;
0313                 lmfdata[index].x = dd.x;
0314                 lmfdata[index].y = dd.y;
0315                 lmfdata[index].z = dd.z;
0316 
0317                 std::list<int>::iterator logicIdIsThere = std::find(mtqLogicIds.begin(), mtqLogicIds.end(), logic_id);
0318                 if (logicIdIsThere != mtqLogicIds.end()) {
0319                   lmfmtq[index].fit_method = mtqConf.getFitMethod(logic_id);
0320                   lmfmtq[index].mtq_ampl = mtqConf.getMTQAmplification(logic_id);
0321                   lmfmtq[index].mtq_time = mtqConf.getMTQTime(logic_id);
0322                   lmfmtq[index].mtq_rise = mtqConf.getMTQRise(logic_id);
0323                   lmfmtq[index].mtq_fwhm = mtqConf.getMTQFWHM(logic_id);
0324                   lmfmtq[index].mtq_fw20 = mtqConf.getMTQFW20(logic_id);
0325                   lmfmtq[index].mtq_fw80 = mtqConf.getMTQFW80(logic_id);
0326                   lmfmtq[index].mtq_sliding = mtqConf.getMTQSliding(logic_id);
0327                 }
0328 
0329                 logicIdIsThere = std::find(laserConfigLogicIds.begin(), laserConfigLogicIds.end(), logic_id);
0330                 if (logicIdIsThere != laserConfigLogicIds.end()) {
0331                   lmflasConf[index].wavelength = laserConfig.getWavelength(logic_id);
0332                   lmflasConf[index].vfe_gain = laserConfig.getVFEGain(logic_id);
0333                   lmflasConf[index].pn_gain = laserConfig.getPNGain(logic_id);
0334                   lmflasConf[index].lsr_power = laserConfig.getLSRPower(logic_id);
0335                   lmflasConf[index].lsr_attenuator = laserConfig.getLSRAttenuator(logic_id);
0336                   lmflasConf[index].lsr_current = laserConfig.getLSRCurrent(logic_id);
0337                   lmflasConf[index].lsr_delay_1 = laserConfig.getLSRDelay1(logic_id);
0338                   lmflasConf[index].lsr_delay_2 = laserConfig.getLSRDelay1(logic_id);
0339                 }
0340                 xcount++;
0341               }
0342               //
0343               li++;
0344               count++;
0345             }
0346             std::cout << "   RunDat: " << std::setw(4) << count << " (" << std::setw(4) << xcount << ")";
0347             count = 0;
0348             xcount = 0;
0349             /*** pnPrim Dat ***/
0350             logic_ids.clear();
0351             logic_ids = pnPrim.getLogicIds();
0352             li = logic_ids.begin();
0353             le = logic_ids.end();
0354             while (li != le) {
0355               int logic_id = *li;
0356               std::vector<EcalLogicID> xtals = econn->getEcalLogicIDForLMPN(logic_id);
0357               for (unsigned int j = 0; j < xtals.size(); j++) {
0358                 int xtal_id = xtals[j].getLogicID();
0359                 int detId = detids[xtal_id];
0360                 detIdData dd = getCoords(detId);
0361                 int index = dd.hashedIndex;
0362                 lmfdata[index].pnMean = pnPrim.getMean(logic_id);
0363                 lmfdata[index].pnRMS = pnPrim.getRMS(logic_id);
0364                 lmfdata[index].pnM3 = pnPrim.getM3(logic_id);
0365                 lmfdata[index].pnAoverpnBMean = pnPrim.getPNAoverBMean(logic_id);
0366                 lmfdata[index].pnAoverpnBRMS = pnPrim.getPNAoverBRMS(logic_id);
0367                 lmfdata[index].pnAoverpnBM3 = pnPrim.getPNAoverBM3(logic_id);
0368                 lmfdata[index].pnFlag = pnPrim.getFlag(logic_id);
0369                 lmfdata[index].detId = detids[xtal_id];
0370                 lmfdata[index].logic_id = xtal_id;
0371                 lmfdata[index].eta = dd.eta;
0372                 lmfdata[index].phi = dd.phi;
0373                 lmfdata[index].x = dd.x;
0374                 lmfdata[index].y = dd.y;
0375                 lmfdata[index].z = dd.z;
0376                 xcount++;
0377               }
0378               li++;
0379               count++;
0380             }
0381             std::cout << "   PnDat: " << std::setw(4) << count << " (" << std::setw(4) << xcount << ")";
0382             count = 0;
0383             xcount = 0;
0384             /*** Prm Dat ***/
0385             logic_ids.clear();
0386             logic_ids = prim.getLogicIds();
0387             li = logic_ids.begin();
0388             le = logic_ids.end();
0389             while (li != le) {
0390               int logic_id = *li;
0391               int detId = detids[logic_id];
0392               detIdData dd = getCoords(detId);
0393               int index = dd.hashedIndex;
0394               lmfdata[index].detId = detId;
0395               lmfdata[index].logic_id = logic_id;
0396               lmfdata[index].eta = dd.eta;
0397               lmfdata[index].phi = dd.phi;
0398               lmfdata[index].x = dd.x;
0399               lmfdata[index].y = dd.y;
0400               lmfdata[index].z = dd.z;
0401               //
0402               lmfdata[index].Mean = prim.getMean(logic_id);
0403               lmfdata[index].RMS = prim.getRMS(logic_id);
0404               lmfdata[index].M3 = prim.getM3(logic_id);
0405               lmfdata[index].APDoverPNAMean = prim.getAPDoverAMean(logic_id);
0406               lmfdata[index].APDoverPNARMS = prim.getAPDoverARMS(logic_id);
0407               lmfdata[index].APDoverPNAM3 = prim.getAPDoverAM3(logic_id);
0408               lmfdata[index].APDoverPNBMean = prim.getAPDoverBMean(logic_id);
0409               lmfdata[index].APDoverPNBRMS = prim.getAPDoverBRMS(logic_id);
0410               lmfdata[index].APDoverPNBM3 = prim.getAPDoverBM3(logic_id);
0411               lmfdata[index].APDoverPNMean = prim.getAPDoverPnMean(logic_id);
0412               lmfdata[index].APDoverPNRMS = prim.getAPDoverPnRMS(logic_id);
0413               lmfdata[index].APDoverPNM3 = prim.getAPDoverPnM3(logic_id);
0414               lmfdata[index].alpha = prim.getAlpha(logic_id);
0415               lmfdata[index].beta = prim.getBeta(logic_id);
0416               lmfdata[index].flag = prim.getFlag(logic_id);
0417               li++;
0418               count++;
0419               xcount++;
0420             }
0421             std::cout << "   PrimDat: " << std::setw(4) << count << " (" << std::setw(4) << xcount << ")";
0422             // fill the tree
0423             xcount = 0;
0424             for (unsigned int i = 0; i < lmfdata.size(); i++) {
0425               if ((userCheck(lmfdata[i])) && (lmfdata[i].Mean > 0)) {
0426                 detId = lmfdata[i].detId;
0427                 logic_id = lmfdata[i].logic_id;
0428                 eta = lmfdata[i].eta;
0429                 phi = lmfdata[i].phi;
0430                 x = lmfdata[i].x;
0431                 y = lmfdata[i].y;
0432                 z = lmfdata[i].z;
0433                 nevts = lmfdata[i].nevts;
0434                 quality = lmfdata[i].quality;
0435                 pnMean = lmfdata[i].pnMean;
0436                 pnRMS = lmfdata[i].pnRMS;
0437                 pnM3 = lmfdata[i].pnM3;
0438                 pnAoverBMean = lmfdata[i].pnAoverpnBMean;
0439                 pnAoverBRMS = lmfdata[i].pnAoverpnBRMS;
0440                 pnAoverBM3 = lmfdata[i].pnAoverpnBM3;
0441                 pnFlag = lmfdata[i].pnFlag;
0442                 Mean = lmfdata[i].Mean;
0443                 RMS = lmfdata[i].RMS;
0444                 M3 = lmfdata[i].M3;
0445                 APDoverPNAMean = lmfdata[i].APDoverPNAMean;
0446                 APDoverPNARMS = lmfdata[i].APDoverPNARMS;
0447                 APDoverPNAM3 = lmfdata[i].APDoverPNAM3;
0448                 APDoverPNBMean = lmfdata[i].APDoverPNBMean;
0449                 APDoverPNBRMS = lmfdata[i].APDoverPNBRMS;
0450                 APDoverPNBM3 = lmfdata[i].APDoverPNBM3;
0451                 APDoverPNMean = lmfdata[i].APDoverPNMean;
0452                 APDoverPNRMS = lmfdata[i].APDoverPNRMS;
0453                 APDoverPNM3 = lmfdata[i].APDoverPNM3;
0454                 alpha = lmfdata[i].alpha;
0455                 beta = lmfdata[i].beta;
0456                 flag = lmfdata[i].flag;
0457                 // matacq variables
0458                 fit_method = lmfmtq[i].fit_method;
0459                 mtq_ampl = lmfmtq[i].mtq_ampl;
0460                 mtq_time = lmfmtq[i].mtq_time;
0461                 mtq_rise = lmfmtq[i].mtq_rise;
0462                 mtq_fwhm = lmfmtq[i].mtq_fwhm;
0463                 mtq_fw20 = lmfmtq[i].mtq_fw20;
0464                 mtq_fw80 = lmfmtq[i].mtq_fw80;
0465                 mtq_sliding = lmfmtq[i].mtq_sliding;
0466                 // laser conf variables
0467                 wavelength = lmflasConf[i].wavelength;
0468                 vfe_gain = lmflasConf[i].vfe_gain;
0469                 pn_gain = lmflasConf[i].pn_gain;
0470                 lsr_power = lmflasConf[i].lsr_power;
0471                 lsr_attenuator = lmflasConf[i].lsr_attenuator;
0472                 lsr_current = lmflasConf[i].lsr_current;
0473                 lsr_delay_1 = lmflasConf[i].lsr_delay_1;
0474                 lsr_delay_2 = lmflasConf[i].lsr_delay_2;
0475 
0476                 tree->Fill();
0477                 xcount++;
0478               }
0479             }
0480             std::cout << " Tree: " << xcount;
0481           }
0482           ri++;
0483         }
0484         si++;
0485       }
0486       std::cout << std::endl;
0487     }
0488   }
0489 
0490 private:
0491   CondDBApp() = delete;  // hidden default constructor
0492   void init();
0493   TTree *tree;
0494   TFile *tfile;
0495   EcalCondDBInterface *econn;
0496   run_t run_min;
0497   run_t run_max;
0498 
0499   int run_num;
0500   uint64_t seqStart;
0501   uint64_t seqStop;
0502   int seqNum;
0503   int lmr;
0504   char color[35];
0505   int nevts;
0506   uint64_t subrunstart;
0507   uint64_t subrunstop;
0508   int detId;
0509   int logic_id;
0510   int eta;
0511   int phi;
0512   int x;
0513   int y;
0514   int z;
0515   int quality;
0516   float pnMean;
0517   float pnRMS;
0518   float pnM3;
0519   float pnAoverBMean;
0520   float pnAoverBRMS;
0521   float pnAoverBM3;
0522   float pnFlag;
0523   float Mean;
0524   float RMS;
0525   float M3;
0526   float APDoverPNAMean;
0527   float APDoverPNARMS;
0528   float APDoverPNAM3;
0529   float APDoverPNBMean;
0530   float APDoverPNBRMS;
0531   float APDoverPNBM3;
0532   float APDoverPNMean;
0533   float APDoverPNRMS;
0534   float APDoverPNM3;
0535   float alpha;
0536   float beta;
0537   int flag;
0538   // matacq
0539   int fit_method;
0540   float mtq_ampl;
0541   float mtq_time;
0542   float mtq_rise;
0543   float mtq_fwhm;
0544   float mtq_fw20;
0545   float mtq_fw80;
0546   float mtq_sliding;
0547   // laser config
0548   int wavelength;
0549   float vfe_gain;
0550   float pn_gain;
0551   float lsr_power;
0552   float lsr_attenuator;
0553   float lsr_current;
0554   float lsr_delay_1;
0555   float lsr_delay_2;
0556 };
0557 
0558 void CondDBApp::init() {
0559   std::stringstream title;
0560   std::stringstream fname;
0561   title << "Dump of Laser data from online DB for run";
0562   fname << "DumpLaserDB";
0563   if (run_max <= 0) {
0564     title << " " << run_min;
0565     fname << "-" << run_min;
0566     run_max = run_min;
0567   } else {
0568     title << "s " << run_min << " - " << run_max;
0569     fname << "-" << run_min << "-" << run_max;
0570   }
0571   fname << ".root";
0572   std::cout << "Building tree " << title.str() << " on file " << fname.str() << std::endl;
0573   tfile = new TFile(fname.str().c_str(), "RECREATE", title.str().c_str());
0574   tree = new TTree("LDB", title.str().c_str());
0575   tree->Branch("run", &run_num, "run/I");
0576   tree->Branch("seqStart", &seqStart, "seqStart/l");
0577   tree->Branch("seqStop", &seqStop, "seqStop/l");
0578   tree->Branch("seqNum", &seqNum, "seqNum/I");
0579   tree->Branch("lmr", &lmr, "lmr/I");
0580   tree->Branch("nevts", &nevts, "nevts/I");
0581   tree->Branch("color", &color, "color/C");
0582   tree->Branch("subrunstart", &subrunstart, "subrunstart/l");
0583   tree->Branch("subrunstop", &subrunstart, "subrunstop/l");
0584   tree->Branch("detId", &detId, "detId/I");
0585   tree->Branch("logic_id", &logic_id, "logic_id/I");
0586   tree->Branch("eta", &eta, "eta/I");
0587   tree->Branch("phi", &phi, "phi/I");
0588   tree->Branch("x", &x, "x/I");
0589   tree->Branch("y", &y, "y/I");
0590   tree->Branch("z", &z, "z/I");
0591   tree->Branch("quality", &quality, "quality/I");
0592   tree->Branch("pnMean", &pnMean, "pnMean/F");
0593   tree->Branch("pnRMS", &pnRMS, "pnRMS/F");
0594   tree->Branch("pnM3", &pnM3, "pnM3/F");
0595   tree->Branch("pnAoverBMean", &pnAoverBMean, "pnAoverBMean/F");
0596   tree->Branch("pnAoverBRMS", &pnAoverBRMS, "pnAoverBRMS/F");
0597   tree->Branch("pnAoverBM3", &pnAoverBM3, "pnAoverBM3/F");
0598   tree->Branch("pnFlag", &pnFlag, "pnFlag/F");
0599   tree->Branch("Mean", &pnMean, "Mean/F");
0600   tree->Branch("RMS", &pnRMS, "RMS/F");
0601   tree->Branch("M3", &pnM3, "M3/F");
0602   tree->Branch("APDoverPnAMean", &APDoverPNAMean, "APDoverPNAMean/F");
0603   tree->Branch("APDoverPnARMS", &APDoverPNARMS, "APDoverPNARMS/F");
0604   tree->Branch("APDoverPnAM3", &APDoverPNAM3, "APDoverPNAM3/F");
0605   tree->Branch("APDoverPnBMean", &APDoverPNBMean, "APDoverPNBMean/F");
0606   tree->Branch("APDoverPnBRMS", &APDoverPNBRMS, "APDoverPNBRMS/F");
0607   tree->Branch("APDoverPnBM3", &APDoverPNBM3, "APDoverPNBM3/F");
0608   tree->Branch("APDoverPnMean", &APDoverPNMean, "APDoverPNMean/F");
0609   tree->Branch("APDoverPnRMS", &APDoverPNRMS, "APDoverPNRMS/F");
0610   tree->Branch("APDoverPnM3", &APDoverPNM3, "APDoverPNM3/F");
0611   tree->Branch("alpha", &alpha, "alpha/F");
0612   tree->Branch("beta", &beta, "beta/F");
0613   tree->Branch("flag", &flag, "flag/I");
0614 
0615   tree->Branch("fit_method", &fit_method, "fit_method/I");
0616   tree->Branch("mtq_ampl", &mtq_ampl, "mtq_ampl/F");
0617   tree->Branch("mtq_time", &mtq_time, "mtq_time/F");
0618   tree->Branch("mtq_rise", &mtq_rise, "mtq_rise/F");
0619   tree->Branch("mtq_fwhm", &mtq_fwhm, "mtq_fwhm/F");
0620   tree->Branch("mtq_fw20", &mtq_fw20, "mtq_fw20/F");
0621   tree->Branch("mtq_fw80", &mtq_fw80, "mtq_fw80/F");
0622   tree->Branch("mtq_sliding", &mtq_sliding, "mtq_sliding/F");
0623 
0624   tree->Branch("wavelength", &wavelength, "wavelength/I");
0625   tree->Branch("vfe_gain", &vfe_gain, "vfe_gain/F");
0626   tree->Branch("pn_gain", &pn_gain, "pn_gain/F");
0627   tree->Branch("lsr_power", &lsr_power, "lsr_power/F");
0628   tree->Branch("lsr_attenuator", &lsr_attenuator, "lsr_attenuator/F");
0629   tree->Branch("lsr_current", &lsr_current, "lsr_current/F");
0630   tree->Branch("lsr_delay_1", &lsr_delay_1, "lsr_delay_1/F");
0631   tree->Branch("lsr_delay_2", &lsr_delay_2, "lsr_delay_2/F");
0632 
0633   std::cout << "Tree created" << std::endl;
0634 }
0635 
0636 void help() {
0637   std::cout << "DumpLaserDB\n";
0638   std::cout << " Reads laser DB and build a root ntuple.\n";
0639   std::cout << std::endl;
0640   std::cout << " Loops on runs between run_min and (if present) run_max.\n";
0641   std::cout << " For each run loops on each Laser Sequence. Each sequence\n";
0642   std::cout << " contains up to 92 LMR.\n";
0643   std::cout << " For each region gets:\n";
0644   std::cout << "   a) run data, such as the number of events acquired\n";
0645   std::cout << "   b) pn primitives\n";
0646   std::cout << "   c) laser primitives\n";
0647   std::cout << std::endl;
0648   std::cout << " During execution it shows the progress on Sequences/LMR's\n";
0649   std::cout << " For each LMR shows its size (crystals belonging to it).\n";
0650   std::cout << " For each group of data shows the number of rows found\n";
0651   std::cout << " in the DB and the no. of channels (in parenthesis).\n";
0652   std::cout << " Usually there is 1 row per RunDat (common to all channels\n";
0653   std::cout << " in the LMR) and as many channels as crystals belonging to\n";
0654   std::cout << " the LMR. There are between between 8 and 10 rows per PnDat\n";
0655   std::cout << " grouped in such a way that the number of channels must be\n";
0656   std::cout << " equal to the crystals in the supermodule. Primitives are \n";
0657   std::cout << " less or equal to the number of channels in the LMR.\n";
0658   std::cout << " Finally, the number of rows in the Tree is shown\n\n";
0659 }
0660 
0661 int main(int argc, char *argv[]) {
0662   namespace po = boost::program_options;
0663 
0664   int run1 = 0;
0665   int run2 = 0;
0666 
0667   std::string sid = "CMS_ORCOFF_PROD";
0668   std::string user = "CMS_ECAL_R";
0669   std::string pass = "";
0670 
0671   std::string confFile = "";
0672 
0673   try {
0674     // options definition
0675     po::options_description desc("Allowed options");
0676     desc.add_options()("help", "shows help message")(
0677         "rmin", po::value<int>()->default_value(0), "set minimum run number to analyze (mandatory)")(
0678         "rmax", po::value<int>()->default_value(0), "set maximum run number to analyze (optional)")(
0679         "sid", po::value<std::string>()->default_value(sid), "Oracle System ID (SID) (defaulted)")(
0680         "user", po::value<std::string>()->default_value(user), "SID user (defaulted)")(
0681         "pass", po::value<std::string>(), "password (mandatory)")(
0682         "file", po::value<std::string>(), "configuration file name (optional)");
0683 
0684     po::positional_options_description p;
0685     p.add("sid", 1);
0686     p.add("user", 1);
0687     p.add("pass", 1);
0688     p.add("rmin", 1);
0689     p.add("rmax", 1);
0690 
0691     // parsing and decoding options
0692     po::variables_map vm;
0693     //    po::store(po::parse_command_line(argc, argv, desc), vm);
0694     po::store(po::command_line_parser(argc, argv).options(desc).positional(p).run(), vm);
0695 
0696     if (vm.count("help")) {
0697       std::cout << desc << std::endl;
0698       help();
0699       return 1;
0700     }
0701 
0702     if (vm.count("rmin") && (vm["rmin"].as<int>() > 0)) {
0703       run1 = vm["rmin"].as<int>();
0704       run2 = vm["rmax"].as<int>();
0705     } else {
0706       std::cout << desc << std::endl;
0707       return 1;
0708     }
0709     if (vm.count("rmax")) {
0710       run2 = vm["rmax"].as<int>();
0711     }
0712     if (vm.count("sid")) {
0713       sid = vm["sid"].as<std::string>();
0714     }
0715     if (vm.count("user")) {
0716       user = vm["user"].as<std::string>();
0717     }
0718     if (vm.count("pass")) {
0719       pass = vm["pass"].as<std::string>();
0720     } else {
0721       std::cout << desc << std::endl;
0722       return 1;
0723     }
0724     if (vm.count("file")) {
0725       confFile = vm["file"].as<std::string>();
0726     }
0727     po::notify(vm);
0728   } catch (std::exception &e) {
0729     std::cout << e.what() << std::endl;
0730     return 1;
0731   }
0732 
0733   try {
0734     CondDBApp app(sid, user, pass, run1, run2);
0735     app.doRun();
0736   } catch (std::exception &e) {
0737     std::cout << "ERROR:  " << e.what() << std::endl;
0738   } catch (...) {
0739     std::cout << "Unknown error caught" << std::endl;
0740   }
0741 
0742   std::cout << "All Done." << std::endl;
0743 
0744   return 0;
0745 }