DTCalibration/plugins/DTTTrigCorrectionFirst.cc

0001
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \author S. Maselli - INFN Torino
0006  */
0007
0008 #include "CalibMuon/DTCalibration/plugins/DTTTrigCorrectionFirst.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0015 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0016 #include "Geometry/DTGeometry/interface/DTSuperLayer.h"
0017
0018 #include "CondFormats/DTObjects/interface/DTTtrig.h"
0019 #include "CondFormats/DataRecord/interface/DTTtrigRcd.h"
0020 #include "CalibMuon/DTCalibration/interface/DTCalibDBUtils.h"
0021
0022 #include <iostream>
0023 #include <fstream>
0024
0025 using namespace edm;
0026 using namespace std;
0027
0028 DTTTrigCorrectionFirst::DTTTrigCorrectionFirst(const ParameterSet& pset)
0029     : dtGeomToken_(esConsumes<edm::Transition::BeginRun>()),
0030       ttrigToken_(
0031           esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", pset.getUntrackedParameter<string>("dbLabel")))) {
0032   debug = pset.getUntrackedParameter<bool>("debug", false);
0033   ttrigMin = pset.getUntrackedParameter<double>("ttrigMin", 0);
0034   ttrigMax = pset.getUntrackedParameter<double>("ttrigMax", 5000);
0035   rmsLimit = pset.getUntrackedParameter<double>("rmsLimit", 5.);
0036 }
0037
0038 DTTTrigCorrectionFirst::~DTTTrigCorrectionFirst() {}
0039
0040 void DTTTrigCorrectionFirst::beginRun(const edm::Run& run, const edm::EventSetup& setup) {
0041   ESHandle<DTTtrig> tTrig;
0042   tTrig = setup.getHandle(ttrigToken_);
0043   tTrigMap = &setup.getData(ttrigToken_);
0044   cout << "[DTTTrigCorrection]: TTrig version: " << tTrig->version() << endl;
0045   muonGeom = setup.getHandle(dtGeomToken_);
0046 }
0047
0048 void DTTTrigCorrectionFirst::endJob() {
0049   // Create the object to be written to DB
0050   DTTtrig tTrigNewMap;
0051   //Get the superlayers list
0052   vector<const DTSuperLayer*> dtSupLylist = muonGeom->superLayers();
0053
0054   //Loop on all superlayers to compute the mean
0055   double average = 0.;
0056   double average2 = 0.;
0057   double rms = 0.;
0058   double averageSigma = 0.;
0059   double counter = 0.;
0060   double averagekfactor = 0;
0061   float kfactor = 0;
0062
0063   for (auto sl = dtSupLylist.begin(); sl != dtSupLylist.end(); sl++) {
0064     float ttrigMean = 0;
0065     float ttrigSigma = 0;
0066     tTrigMap->get((*sl)->id(), ttrigMean, ttrigSigma, kfactor, DTTimeUnits::ns);
0067     if (ttrigMean < ttrigMax && ttrigMean > ttrigMin) {
0068       average += ttrigMean;
0069       averageSigma += ttrigSigma;
0070       if (ttrigMean > 0)
0071         counter += 1.;
0072     }
0073   }  //End of loop on superlayers
0074
0075   average = average / counter;
0076   averageSigma = averageSigma / counter;
0077   averagekfactor = kfactor;
0078
0079   //  cout << " average counter "<< average << " "<< counter <<endl;
0080
0081   for (auto sl = dtSupLylist.begin(); sl != dtSupLylist.end(); sl++) {
0082     float ttrigMean = 0;
0083     float ttrigSigma = 0;
0084     float kfactor = 0;
0085     tTrigMap->get((*sl)->id(), ttrigMean, ttrigSigma, kfactor, DTTimeUnits::ns);
0086     if (ttrigMean < ttrigMax && ttrigMean > ttrigMin) {
0087       average2 += (ttrigMean - average) * (ttrigMean - average);
0088     }
0089   }  //End of loop on superlayers
0090
0091   rms = average2 / (counter - 1);
0092   rms = sqrt(rms);
0093   cout << "average averageSigma counter rms " << average << " " << averageSigma << " " << counter << " " << rms << endl;
0094
0095   for (auto sl = dtSupLylist.begin(); sl != dtSupLylist.end(); sl++) {
0096     //Compute new ttrig
0097     double newTTrigMean = 0;
0098     double newTTrigSigma = 0;
0099     double newkfactor = 0;
0100     float tempttrigMean = 0;
0101     float tempttrigSigma = 0;
0102     float tempkfactor = 0;
0103     float ttrigMean = 0;
0104     float ttrigSigma = 0;
0105     float kfactor = 0;
0106     tTrigMap->get((*sl)->id(), ttrigMean, ttrigSigma, kfactor, DTTimeUnits::ns);
0107     int chamber = (*sl)->id().chamberId().station();
0108
0109     //check if ttrigMean is similar to the mean
0110     if (abs(ttrigMean - average) < rmsLimit * rms) {
0111       newTTrigMean = ttrigMean;
0112       newTTrigSigma = ttrigSigma;
0113       newkfactor = kfactor;
0114     } else {
0115       // do not consider if ttrig == 0
0116       if (ttrigMean > 0) {
0117         //cout << "ttrig chamber " << ttrigMean <<" "<<chamber<<endl;
0118         if (((*sl)->id().superlayer()) == 1) {
0119           //cout << " superlayer " << ((*sl)->id().superlayer()) <<endl;
0120           DTSuperLayerId slId((*sl)->id().chamberId(), 3);
0121           tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0122           if (abs(tempttrigMean - average) < rmsLimit * rms) {
0123             newTTrigMean = tempttrigMean;
0124             newTTrigSigma = tempttrigSigma;
0125             newkfactor = tempkfactor;
0126             cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0127                  << " -> takes value of sl 3 " << newTTrigMean << endl;
0128           } else if (chamber == 4) {
0129             cout << "No correction possible within same chamber (sl1) " << endl;
0130             newTTrigMean = average;
0131             newTTrigSigma = averageSigma;
0132             newkfactor = averagekfactor;
0133             cout << "####### Bad SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0134           } else if (chamber != 4) {
0135             DTSuperLayerId slId((*sl)->id().chamberId(), 2);
0136             tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0137             if (abs(tempttrigMean - average) < rmsLimit * rms) {
0138               newTTrigMean = tempttrigMean;
0139               newTTrigSigma = tempttrigSigma;
0140               newkfactor = tempkfactor;
0141               cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0142                    << " -> takes value of sl 2 " << newTTrigMean << endl;
0143             } else {
0144               cout << "No correction possible within same chamber (sl1) " << endl;
0145               newTTrigMean = average;
0146               newTTrigSigma = averageSigma;
0147               newkfactor = averagekfactor;
0148               cout << "####### Bad SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0149             }
0150           }
0151         } else if (((*sl)->id().superlayer()) == 2) {
0152           //cout << " superlayer " << ((*sl)->id().superlayer()) <<endl;
0153           DTSuperLayerId slId((*sl)->id().chamberId(), 1);
0154           tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0155           if (abs(tempttrigMean - average) < rmsLimit * rms) {
0156             newTTrigMean = tempttrigMean;
0157             newTTrigSigma = tempttrigSigma;
0158             newkfactor = tempkfactor;
0159             cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0160                  << " -> takes value of sl 1 " << newTTrigMean << endl;
0161           } else {
0162             DTSuperLayerId slId((*sl)->id().chamberId(), 3);
0163             tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0164             if (abs(tempttrigMean - average) < rmsLimit * rms) {
0165               newTTrigMean = tempttrigMean;
0166               newTTrigSigma = tempttrigSigma;
0167               newkfactor = tempkfactor;
0168               cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0169                    << "-> takes value of sl 3 " << newTTrigMean << endl;
0170             } else {
0171               cout << "No correction possible within same chamber (sl2)  " << endl;
0172               newTTrigMean = average;
0173               newTTrigSigma = averageSigma;
0174               newkfactor = averagekfactor;
0175               cout << "####### Bad SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0176             }
0177           }
0178
0179         } else if (((*sl)->id().superlayer()) == 3) {
0180           //cout << " superlayer " << ((*sl)->id().superlayer()) <<endl;
0181           DTSuperLayerId slId((*sl)->id().chamberId(), 1);
0182           tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0183           if (abs(tempttrigMean - average) < rmsLimit * rms) {
0184             newTTrigMean = tempttrigMean;
0185             newTTrigSigma = tempttrigSigma;
0186             newkfactor = tempkfactor;
0187             cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0188                  << " -> takes value of sl 1 " << newTTrigMean << endl;
0189           } else if (chamber == 4) {
0190             cout << "No correction possible within same chamber (sl3)" << endl;
0191             newTTrigMean = average;
0192             newTTrigSigma = averageSigma;
0193             newkfactor = averagekfactor;
0194             cout << "####### Bad SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0195           } else if (chamber != 4) {
0196             DTSuperLayerId slId((*sl)->id().chamberId(), 2);
0197             tTrigMap->get(slId, tempttrigMean, tempttrigSigma, tempkfactor, DTTimeUnits::ns);
0198             if (abs(tempttrigMean - average) < rmsLimit * rms) {
0199               newTTrigMean = tempttrigMean;
0200               newTTrigSigma = tempttrigSigma;
0201               newkfactor = tempkfactor;
0202               cout << "Chamber " << chamber << " sl " << ((*sl)->id().superlayer()) << "has ttrig " << ttrigMean
0203                    << " -> takes value of sl 2 " << newTTrigMean << endl;
0204             } else {
0205               cout << "No correction possible within same chamber (sl3) " << endl;
0206               newTTrigMean = average;
0207               newTTrigSigma = averageSigma;
0208               newkfactor = averagekfactor;
0209               cout << "####### Bad SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0210             }
0211           }
0212         }
0213       } else {
0214         //  cout << "SL not present " << ttrigMean << " " <<((*sl)->id().superlayer()) <<endl;
0215         newTTrigMean = average;
0216         newTTrigSigma = averageSigma;
0217         newkfactor = averagekfactor;
0218         cout << "####### NotPresent SL: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0219       }
0220     }
0221
0222     //Store new ttrig in the new map
0223     tTrigNewMap.set((*sl)->id(), newTTrigMean, newTTrigSigma, newkfactor, DTTimeUnits::ns);
0224     if (debug) {
0225       cout << "New tTrig: " << (*sl)->id() << " from " << ttrigMean << " to " << newTTrigMean << endl;
0226       cout << "New tTrigSigma: " << (*sl)->id() << " from " << ttrigSigma << " to " << newTTrigSigma << endl;
0227     }
0228   }  //End of loop on superlayers
0229
0230   //Write object to DB
0231   cout << "[DTTTrigCorrection]: Writing ttrig object to DB!" << endl;
0232   string record = "DTTtrigRcd";
0233   DTCalibDBUtils::writeToDB<DTTtrig>(record, tTrigNewMap);
0234 }