00001 #include "HcalHitReconstructor.h"
00002 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
00003 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
00004 #include "DataFormats/Common/interface/EDCollection.h"
00005 #include "DataFormats/Common/interface/Handle.h"
00006 #include "FWCore/Framework/interface/Selector.h"
00007 #include "FWCore/Framework/interface/ESHandle.h"
00008 #include "FWCore/Framework/interface/EventSetup.h"
00009 #include "CalibFormats/HcalObjects/interface/HcalCoderDb.h"
00010 #include "CalibFormats/HcalObjects/interface/HcalCalibrations.h"
00011 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
00012 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
00013 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputer.h"
00014 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputerRcd.h"
00015 #include "CalibCalorimetry/HcalAlgos/interface/HcalDbASCIIIO.h"
00016
00017 #include <iostream>
00018 #include <fstream>
00019
00020
00021
00022
00023 HcalHitReconstructor::HcalHitReconstructor(edm::ParameterSet const& conf):
00024 reco_(conf.getParameter<bool>("correctForTimeslew"),
00025 conf.getParameter<bool>("correctForPhaseContainment"),
00026 conf.getParameter<double>("correctionPhaseNS")),
00027 det_(DetId::Hcal),
00028 inputLabel_(conf.getParameter<edm::InputTag>("digiLabel")),
00029 correctTiming_(conf.getParameter<bool>("correctTiming")),
00030 setNoiseFlags_(conf.getParameter<bool>("setNoiseFlags")),
00031 setHSCPFlags_(conf.getParameter<bool>("setHSCPFlags")),
00032 setSaturationFlags_(conf.getParameter<bool>("setSaturationFlags")),
00033 setTimingTrustFlags_(conf.getParameter<bool>("setTimingTrustFlags")),
00034 setPulseShapeFlags_(conf.getParameter<bool>("setPulseShapeFlags")),
00035 dropZSmarkedPassed_(conf.getParameter<bool>("dropZSmarkedPassed")),
00036 firstAuxTS_(conf.getParameter<int>("firstAuxTS")),
00037 firstSample_(conf.getParameter<int>("firstSample")),
00038 samplesToAdd_(conf.getParameter<int>("samplesToAdd")),
00039 tsFromDB_(conf.getParameter<bool>("tsFromDB")),
00040 useLeakCorrection_( conf.getParameter<bool>("useLeakCorrection"))
00041 {
00042
00043 std::string subd=conf.getParameter<std::string>("Subdetector");
00044
00045
00046
00047
00048
00049 recoParamsFromDB_ = conf.getParameter<bool>("recoParamsFromDB");
00050
00051
00052
00053
00054 hbheFlagSetter_ = 0;
00055 hbheHSCPFlagSetter_ = 0;
00056 hbhePulseShapeFlagSetter_ = 0;
00057 hbheTimingShapedFlagSetter_ = 0;
00058 hfdigibit_ = 0;
00059
00060 hfS9S1_ = 0;
00061 hfS8S1_ = 0;
00062 hfPET_ = 0;
00063 saturationFlagSetter_ = 0;
00064 HFTimingTrustFlagSetter_ = 0;
00065 digiTimeFromDB_ = false;
00066
00067 if (setSaturationFlags_)
00068 {
00069 const edm::ParameterSet& pssat = conf.getParameter<edm::ParameterSet>("saturationParameters");
00070 saturationFlagSetter_ = new HcalADCSaturationFlag(pssat.getParameter<int>("maxADCvalue"));
00071 }
00072
00073 if (!strcasecmp(subd.c_str(),"HBHE")) {
00074 subdet_=HcalBarrel;
00075 bool timingShapedCutsFlags = conf.getParameter<bool>("setTimingShapedCutsFlags");
00076 if (timingShapedCutsFlags)
00077 {
00078 const edm::ParameterSet& psTshaped = conf.getParameter<edm::ParameterSet>("timingshapedcutsParameters");
00079 hbheTimingShapedFlagSetter_ = new HBHETimingShapedFlagSetter(psTshaped.getParameter<std::vector<double> >("tfilterEnvelope"),
00080 psTshaped.getParameter<bool>("ignorelowest"),
00081 psTshaped.getParameter<bool>("ignorehighest"),
00082 psTshaped.getParameter<double>("win_offset"),
00083 psTshaped.getParameter<double>("win_gain"));
00084 }
00085
00086 if (setNoiseFlags_)
00087 {
00088 const edm::ParameterSet& psdigi =conf.getParameter<edm::ParameterSet>("flagParameters");
00089 hbheFlagSetter_=new HBHEStatusBitSetter(psdigi.getParameter<double>("nominalPedestal"),
00090 psdigi.getParameter<double>("hitEnergyMinimum"),
00091 psdigi.getParameter<int>("hitMultiplicityThreshold"),
00092 psdigi.getParameter<std::vector<edm::ParameterSet> >("pulseShapeParameterSets")
00093 );
00094 }
00095 if (setHSCPFlags_)
00096 {
00097 const edm::ParameterSet& psHSCP = conf.getParameter<edm::ParameterSet>("hscpParameters");
00098 hbheHSCPFlagSetter_ = new HBHETimeProfileStatusBitSetter(psHSCP.getParameter<double>("r1Min"),
00099 psHSCP.getParameter<double>("r1Max"),
00100 psHSCP.getParameter<double>("r2Min"),
00101 psHSCP.getParameter<double>("r2Max"),
00102 psHSCP.getParameter<double>("fracLeaderMin"),
00103 psHSCP.getParameter<double>("fracLeaderMax"),
00104 psHSCP.getParameter<double>("slopeMin"),
00105 psHSCP.getParameter<double>("slopeMax"),
00106 psHSCP.getParameter<double>("outerMin"),
00107 psHSCP.getParameter<double>("outerMax"),
00108 psHSCP.getParameter<double>("TimingEnergyThreshold"));
00109 }
00110 if (setPulseShapeFlags_)
00111 {
00112 const edm::ParameterSet &psPulseShape = conf.getParameter<edm::ParameterSet>("pulseShapeParameters");
00113 hbhePulseShapeFlagSetter_ = new HBHEPulseShapeFlagSetter(
00114 psPulseShape.getParameter<double>("MinimumChargeThreshold"),
00115 psPulseShape.getParameter<double>("TS4TS5ChargeThreshold"),
00116 psPulseShape.getParameter<unsigned int>("TrianglePeakTS"),
00117 psPulseShape.getParameter<std::vector<double> >("LinearThreshold"),
00118 psPulseShape.getParameter<std::vector<double> >("LinearCut"),
00119 psPulseShape.getParameter<std::vector<double> >("RMS8MaxThreshold"),
00120 psPulseShape.getParameter<std::vector<double> >("RMS8MaxCut"),
00121 psPulseShape.getParameter<std::vector<double> >("LeftSlopeThreshold"),
00122 psPulseShape.getParameter<std::vector<double> >("LeftSlopeCut"),
00123 psPulseShape.getParameter<std::vector<double> >("RightSlopeThreshold"),
00124 psPulseShape.getParameter<std::vector<double> >("RightSlopeCut"),
00125 psPulseShape.getParameter<std::vector<double> >("RightSlopeSmallThreshold"),
00126 psPulseShape.getParameter<std::vector<double> >("RightSlopeSmallCut"),
00127 psPulseShape.getParameter<std::vector<double> >("TS4TS5LowerThreshold"),
00128 psPulseShape.getParameter<std::vector<double> >("TS4TS5LowerCut"),
00129 psPulseShape.getParameter<std::vector<double> >("TS4TS5UpperThreshold"),
00130 psPulseShape.getParameter<std::vector<double> >("TS4TS5UpperCut"),
00131 psPulseShape.getParameter<bool>("UseDualFit"),
00132 psPulseShape.getParameter<bool>("TriangleIgnoreSlow"));
00133 }
00134
00135 produces<HBHERecHitCollection>();
00136 } else if (!strcasecmp(subd.c_str(),"HO")) {
00137 subdet_=HcalOuter;
00138 produces<HORecHitCollection>();
00139 } else if (!strcasecmp(subd.c_str(),"HF")) {
00140 subdet_=HcalForward;
00141 digiTimeFromDB_=conf.getParameter<bool>("digiTimeFromDB");
00142
00143 if (setTimingTrustFlags_) {
00144
00145 const edm::ParameterSet& pstrust = conf.getParameter<edm::ParameterSet>("hfTimingTrustParameters");
00146 HFTimingTrustFlagSetter_=new HFTimingTrustFlag(pstrust.getParameter<int>("hfTimingTrustLevel1"),
00147 pstrust.getParameter<int>("hfTimingTrustLevel2"));
00148 }
00149
00150 if (setNoiseFlags_)
00151 {
00152 const edm::ParameterSet& psdigi =conf.getParameter<edm::ParameterSet>("digistat");
00153 const edm::ParameterSet& psTimeWin =conf.getParameter<edm::ParameterSet>("HFInWindowStat");
00154 hfdigibit_=new HcalHFStatusBitFromDigis(psdigi,psTimeWin);
00155
00156 const edm::ParameterSet& psS9S1 = conf.getParameter<edm::ParameterSet>("S9S1stat");
00157 hfS9S1_ = new HcalHF_S9S1algorithm(psS9S1.getParameter<std::vector<double> >("short_optimumSlope"),
00158 psS9S1.getParameter<std::vector<double> >("shortEnergyParams"),
00159 psS9S1.getParameter<std::vector<double> >("shortETParams"),
00160 psS9S1.getParameter<std::vector<double> >("long_optimumSlope"),
00161 psS9S1.getParameter<std::vector<double> >("longEnergyParams"),
00162 psS9S1.getParameter<std::vector<double> >("longETParams"),
00163 psS9S1.getParameter<int>("HcalAcceptSeverityLevel"),
00164 psS9S1.getParameter<bool>("isS8S1")
00165 );
00166
00167 const edm::ParameterSet& psS8S1 = conf.getParameter<edm::ParameterSet>("S8S1stat");
00168 hfS8S1_ = new HcalHF_S9S1algorithm(psS8S1.getParameter<std::vector<double> >("short_optimumSlope"),
00169 psS8S1.getParameter<std::vector<double> >("shortEnergyParams"),
00170 psS8S1.getParameter<std::vector<double> >("shortETParams"),
00171 psS8S1.getParameter<std::vector<double> >("long_optimumSlope"),
00172 psS8S1.getParameter<std::vector<double> >("longEnergyParams"),
00173 psS8S1.getParameter<std::vector<double> >("longETParams"),
00174 psS8S1.getParameter<int>("HcalAcceptSeverityLevel"),
00175 psS8S1.getParameter<bool>("isS8S1")
00176 );
00177
00178 const edm::ParameterSet& psPET = conf.getParameter<edm::ParameterSet>("PETstat");
00179 hfPET_ = new HcalHF_PETalgorithm(psPET.getParameter<std::vector<double> >("short_R"),
00180 psPET.getParameter<std::vector<double> >("shortEnergyParams"),
00181 psPET.getParameter<std::vector<double> >("shortETParams"),
00182 psPET.getParameter<std::vector<double> >("long_R"),
00183 psPET.getParameter<std::vector<double> >("longEnergyParams"),
00184 psPET.getParameter<std::vector<double> >("longETParams"),
00185 psPET.getParameter<int>("HcalAcceptSeverityLevel"),
00186 psPET.getParameter<std::vector<double> >("short_R_29"),
00187 psPET.getParameter<std::vector<double> >("long_R_29")
00188 );
00189 }
00190 produces<HFRecHitCollection>();
00191 } else if (!strcasecmp(subd.c_str(),"ZDC")) {
00192 det_=DetId::Calo;
00193 subdet_=HcalZDCDetId::SubdetectorId;
00194 produces<ZDCRecHitCollection>();
00195 } else if (!strcasecmp(subd.c_str(),"CALIB")) {
00196 subdet_=HcalOther;
00197 subdetOther_=HcalCalibration;
00198 produces<HcalCalibRecHitCollection>();
00199 } else {
00200 std::cout << "HcalHitReconstructor is not associated with a specific subdetector!" << std::endl;
00201 }
00202
00203 }
00204
00205 HcalHitReconstructor::~HcalHitReconstructor() {
00206 if (hbheFlagSetter_) delete hbheFlagSetter_;
00207 if (hfdigibit_) delete hfdigibit_;
00208 if (hbheHSCPFlagSetter_) delete hbheHSCPFlagSetter_;
00209 if (hbhePulseShapeFlagSetter_) delete hbhePulseShapeFlagSetter_;
00210 if (hfS9S1_) delete hfS9S1_;
00211 if (hfPET_) delete hfPET_;
00212 }
00213
00214 void HcalHitReconstructor::beginRun(edm::Run&r, edm::EventSetup const & es){
00215
00216 if ( tsFromDB_== true || recoParamsFromDB_ == true )
00217 {
00218 edm::ESHandle<HcalRecoParams> p;
00219 es.get<HcalRecoParamsRcd>().get(p);
00220 paramTS = new HcalRecoParams(*p.product());
00221
00222
00223
00224
00225 }
00226
00227 if (digiTimeFromDB_==true)
00228 {
00229 edm::ESHandle<HcalFlagHFDigiTimeParams> p;
00230 es.get<HcalFlagHFDigiTimeParamsRcd>().get(p);
00231 HFDigiTimeParams = new HcalFlagHFDigiTimeParams(*p.product());
00232 }
00233 reco_.beginRun(es);
00234 }
00235
00236 void HcalHitReconstructor::endRun(edm::Run&r, edm::EventSetup const & es){
00237 if (tsFromDB_==true)
00238 {
00239 if (paramTS) delete paramTS;
00240 }
00241 if (digiTimeFromDB_==true)
00242 {
00243 if (HFDigiTimeParams) delete HFDigiTimeParams;
00244 }
00245 reco_.endRun();
00246 }
00247
00248 void HcalHitReconstructor::produce(edm::Event& e, const edm::EventSetup& eventSetup)
00249 {
00250
00251
00252 edm::ESHandle<HcalDbService> conditions;
00253 eventSetup.get<HcalDbRecord>().get(conditions);
00254 const HcalQIEShape* shape = conditions->getHcalShape ();
00255
00256 if(e.isRealData()) reco_.setForData();
00257 if(useLeakCorrection_) reco_.setLeakCorrection();
00258
00259 edm::ESHandle<HcalChannelQuality> p;
00260 eventSetup.get<HcalChannelQualityRcd>().get(p);
00261 HcalChannelQuality* myqual = new HcalChannelQuality(*p.product());
00262
00263 edm::ESHandle<HcalSeverityLevelComputer> mycomputer;
00264 eventSetup.get<HcalSeverityLevelComputerRcd>().get(mycomputer);
00265 const HcalSeverityLevelComputer* mySeverity = mycomputer.product();
00266
00267 if (det_==DetId::Hcal) {
00268
00269
00270 if (subdet_==HcalBarrel || subdet_==HcalEndcap) {
00271 edm::Handle<HBHEDigiCollection> digi;
00272
00273 e.getByLabel(inputLabel_,digi);
00274
00275
00276 std::auto_ptr<HBHERecHitCollection> rec(new HBHERecHitCollection);
00277 rec->reserve(digi->size());
00278
00279 if (setNoiseFlags_) hbheFlagSetter_->Clear();
00280 HBHEDigiCollection::const_iterator i;
00281 std::vector<HBHEDataFrame> HBDigis;
00282 std::vector<int> RecHitIndex;
00283
00284
00285 int favorite_capid = 0;
00286 if (correctTiming_) {
00287 long capid_votes[4] = {0,0,0,0};
00288 for (i=digi->begin(); i!=digi->end(); i++) {
00289 capid_votes[(*i)[0].capid()]++;
00290 }
00291 for (int k = 0; k < 4; k++)
00292 if (capid_votes[k] > capid_votes[favorite_capid])
00293 favorite_capid = k;
00294 }
00295
00296 int first = firstSample_;
00297 int toadd = samplesToAdd_;
00298
00299 for (i=digi->begin(); i!=digi->end(); i++) {
00300 HcalDetId cell = i->id();
00301 DetId detcell=(DetId)cell;
00302
00303 if(tsFromDB_ || recoParamsFromDB_) {
00304 const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00305 firstSample_ = param_ts->firstSample();
00306 samplesToAdd_ = param_ts->samplesToAdd();
00307 if(recoParamsFromDB_) {
00308 bool correctForTimeslew=param_ts->correctForTimeslew();
00309 bool correctForPhaseContainment= param_ts->correctForPhaseContainment();
00310 float phaseNS=param_ts->correctionPhaseNS();
00311 useLeakCorrection_= param_ts->useLeakCorrection();
00312 correctTiming_ = param_ts->correctTiming();
00313 firstAuxTS_ = param_ts->firstAuxTS();
00314 int pileupCleaningID = param_ts->pileupCleaningID();
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341 reco_.setRecoParams(correctForTimeslew,correctForPhaseContainment,useLeakCorrection_,pileupCleaningID,phaseNS);
00342 }
00343 }
00344
00345 first = firstSample_;
00346 toadd = samplesToAdd_;
00347
00348
00349 const HcalChannelStatus* mydigistatus=myqual->getValues(detcell.rawId());
00350 if (mySeverity->dropChannel(mydigistatus->getValue() ) ) continue;
00351 if (dropZSmarkedPassed_)
00352 if (i->zsMarkAndPass()) continue;
00353
00354 const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00355 const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00356 HcalCoderDb coder (*channelCoder, *shape);
00357
00358 rec->push_back(reco_.reconstruct(*i,first,toadd,coder,calibrations));
00359
00360
00361 int auxflag=0;
00362 int fTS = firstAuxTS_;
00363 if (fTS<0) fTS=0;
00364 for (int xx=fTS; xx<fTS+4 && xx<i->size();++xx)
00365 auxflag+=(i->sample(xx).adc())<<(7*(xx-fTS));
00366
00367 auxflag+=((i->sample(fTS).capid())<<28);
00368 (rec->back()).setAux(auxflag);
00369
00370 (rec->back()).setFlags(0);
00371
00372 if (fTS>0)
00373 (rec->back()).setFlagField((i->sample(fTS-1).adc()), HcalCaloFlagLabels::PresampleADC,7);
00374
00375 if (hbheTimingShapedFlagSetter_!=0)
00376 hbheTimingShapedFlagSetter_->SetTimingShapedFlags(rec->back());
00377 if (setNoiseFlags_)
00378 hbheFlagSetter_->SetFlagsFromDigi(rec->back(),*i,coder,calibrations,first,toadd);
00379 if (setPulseShapeFlags_ == true)
00380 hbhePulseShapeFlagSetter_->SetPulseShapeFlags(rec->back(), *i, coder, calibrations);
00381 if (setSaturationFlags_)
00382 saturationFlagSetter_->setSaturationFlag(rec->back(),*i);
00383 if (correctTiming_)
00384 HcalTimingCorrector::Correct(rec->back(), *i, favorite_capid);
00385 if (setHSCPFlags_ && i->id().ietaAbs()<16)
00386 {
00387 double DigiEnergy=0;
00388 for(int j=0; j!=i->size(); DigiEnergy += i->sample(j++).nominal_fC());
00389 if(DigiEnergy > hbheHSCPFlagSetter_->EnergyThreshold())
00390 {
00391 HBDigis.push_back(*i);
00392 RecHitIndex.push_back(rec->size()-1);
00393 }
00394
00395 }
00396 }
00397
00398
00399 if (setNoiseFlags_) hbheFlagSetter_->SetFlagsFromRecHits(*rec);
00400 if (setHSCPFlags_) hbheHSCPFlagSetter_->hbheSetTimeFlagsFromDigi(rec.get(), HBDigis, RecHitIndex);
00401
00402 e.put(rec);
00403
00404
00405 } else if (subdet_==HcalOuter) {
00406 edm::Handle<HODigiCollection> digi;
00407 e.getByLabel(inputLabel_,digi);
00408
00409
00410 std::auto_ptr<HORecHitCollection> rec(new HORecHitCollection);
00411 rec->reserve(digi->size());
00412
00413 HODigiCollection::const_iterator i;
00414
00415
00416 int favorite_capid = 0;
00417 if (correctTiming_) {
00418 long capid_votes[4] = {0,0,0,0};
00419 for (i=digi->begin(); i!=digi->end(); i++) {
00420 capid_votes[(*i)[0].capid()]++;
00421 }
00422 for (int k = 0; k < 4; k++)
00423 if (capid_votes[k] > capid_votes[favorite_capid])
00424 favorite_capid = k;
00425 }
00426
00427 int first = firstSample_;
00428 int toadd = samplesToAdd_;
00429
00430 for (i=digi->begin(); i!=digi->end(); i++) {
00431 HcalDetId cell = i->id();
00432 DetId detcell=(DetId)cell;
00433
00434 if(tsFromDB_ || recoParamsFromDB_) {
00435 const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00436 firstSample_ = param_ts->firstSample();
00437 samplesToAdd_ = param_ts->samplesToAdd();
00438 if(recoParamsFromDB_) {
00439 bool correctForTimeslew=param_ts->correctForTimeslew();
00440 bool correctForPhaseContainment= param_ts->correctForPhaseContainment();
00441 float phaseNS=param_ts->correctionPhaseNS();
00442 useLeakCorrection_= param_ts->useLeakCorrection();
00443 correctTiming_ = param_ts->correctTiming();
00444 firstAuxTS_ = param_ts->firstAuxTS();
00445 int pileupCleaningID = param_ts->pileupCleaningID();
00446 reco_.setRecoParams(correctForTimeslew,correctForPhaseContainment,useLeakCorrection_,pileupCleaningID,phaseNS);
00447 }
00448 }
00449
00450 first = firstSample_;
00451 toadd = samplesToAdd_;
00452
00453
00454 const HcalChannelStatus* mydigistatus=myqual->getValues(detcell.rawId());
00455 if (mySeverity->dropChannel(mydigistatus->getValue() ) ) continue;
00456 if (dropZSmarkedPassed_)
00457 if (i->zsMarkAndPass()) continue;
00458
00459 const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00460 const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00461 HcalCoderDb coder (*channelCoder, *shape);
00462
00463 rec->push_back(reco_.reconstruct(*i,first,toadd,coder,calibrations));
00464
00465
00466 int auxflag=0;
00467 int fTS = firstAuxTS_;
00468 if (fTS<0) fTS=0;
00469 for (int xx=fTS; xx<fTS+4 && xx<i->size();++xx)
00470 auxflag+=(i->sample(xx).adc())<<(7*(xx-fTS));
00471
00472 auxflag+=((i->sample(fTS).capid())<<28);
00473 (rec->back()).setAux(auxflag);
00474
00475 (rec->back()).setFlags(0);
00476
00477 if (fTS>0)
00478 (rec->back()).setFlagField((i->sample(fTS-1).adc()), HcalCaloFlagLabels::PresampleADC,7);
00479
00480 if (setSaturationFlags_)
00481 saturationFlagSetter_->setSaturationFlag(rec->back(),*i);
00482 if (correctTiming_)
00483 HcalTimingCorrector::Correct(rec->back(), *i, favorite_capid);
00484 }
00485
00486 e.put(rec);
00487
00488
00489 } else if (subdet_==HcalForward) {
00490 edm::Handle<HFDigiCollection> digi;
00491 e.getByLabel(inputLabel_,digi);
00492
00493
00495
00496 std::auto_ptr<HFRecHitCollection> rec(new HFRecHitCollection);
00497 rec->reserve(digi->size());
00498
00499 HFDigiCollection::const_iterator i;
00500
00501
00502 int favorite_capid = 0;
00503 if (correctTiming_) {
00504 long capid_votes[4] = {0,0,0,0};
00505 for (i=digi->begin(); i!=digi->end(); i++) {
00506 capid_votes[(*i)[0].capid()]++;
00507 }
00508 for (int k = 0; k < 4; k++)
00509 if (capid_votes[k] > capid_votes[favorite_capid])
00510 favorite_capid = k;
00511 }
00512
00513 int first = firstSample_;
00514 int toadd = samplesToAdd_;
00515
00516 for (i=digi->begin(); i!=digi->end(); i++) {
00517 HcalDetId cell = i->id();
00518 DetId detcell=(DetId)cell;
00519
00520 if(tsFromDB_ || recoParamsFromDB_) {
00521 const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00522 firstSample_ = param_ts->firstSample();
00523 samplesToAdd_ = param_ts->samplesToAdd();
00524 if(recoParamsFromDB_) {
00525 bool correctForTimeslew=param_ts->correctForTimeslew();
00526 bool correctForPhaseContainment= param_ts->correctForPhaseContainment();
00527 float phaseNS=param_ts->correctionPhaseNS();
00528 useLeakCorrection_= param_ts->useLeakCorrection();
00529 correctTiming_ = param_ts->correctTiming();
00530 firstAuxTS_ = param_ts->firstAuxTS();
00531 int pileupCleaningID = param_ts->pileupCleaningID();
00532 reco_.setRecoParams(correctForTimeslew,correctForPhaseContainment,useLeakCorrection_,pileupCleaningID,phaseNS);
00533 }
00534 }
00535
00536 first = firstSample_;
00537 toadd = samplesToAdd_;
00538
00539
00540 const HcalChannelStatus* mydigistatus=myqual->getValues(detcell.rawId());
00541 if (mySeverity->dropChannel(mydigistatus->getValue() ) ) continue;
00542 if (dropZSmarkedPassed_)
00543 if (i->zsMarkAndPass()) continue;
00544
00545 const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00546 const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00547 HcalCoderDb coder (*channelCoder, *shape);
00548
00549
00550 if (digiTimeFromDB_==true && hfdigibit_!=0)
00551 {
00552 const HcalFlagHFDigiTimeParam* hfDTparam = HFDigiTimeParams->getValues(detcell.rawId());
00553 hfdigibit_->resetParamsFromDB(hfDTparam->HFdigiflagFirstSample(),
00554 hfDTparam->HFdigiflagSamplesToAdd(),
00555 hfDTparam->HFdigiflagExpectedPeak(),
00556 hfDTparam->HFdigiflagMinEThreshold(),
00557 hfDTparam->HFdigiflagCoefficients()
00558 );
00559 }
00560
00561
00562 rec->push_back(reco_.reconstruct(*i,first,toadd,coder,calibrations));
00563
00564
00565 int auxflag=0;
00566 int fTS = firstAuxTS_;
00567 if (fTS<0) fTS=0;
00568 for (int xx=fTS; xx<fTS+4 && xx<i->size();++xx)
00569 auxflag+=(i->sample(xx).adc())<<(7*(xx-fTS));
00570
00571 auxflag+=((i->sample(fTS).capid())<<28);
00572 (rec->back()).setAux(auxflag);
00573
00574
00575 (rec->back()).setFlags(0);
00576
00577
00578 if (fTS>0)
00579 (rec->back()).setFlagField((i->sample(fTS-1).adc()), HcalCaloFlagLabels::PresampleADC,7);
00580
00581
00582 if (setNoiseFlags_)
00583 hfdigibit_->hfSetFlagFromDigi(rec->back(),*i,coder,calibrations);
00584 if (setSaturationFlags_)
00585 saturationFlagSetter_->setSaturationFlag(rec->back(),*i);
00586 if (setTimingTrustFlags_)
00587 HFTimingTrustFlagSetter_->setHFTimingTrustFlag(rec->back(),*i);
00588 if (correctTiming_)
00589 HcalTimingCorrector::Correct(rec->back(), *i, favorite_capid);
00590 }
00591
00592
00593
00594
00595 if (setNoiseFlags_)
00596 {
00597
00598
00599
00600 for (HFRecHitCollection::iterator i = rec->begin();i!=rec->end();++i)
00601 {
00602 int depth=i->id().depth();
00603 int ieta=i->id().ieta();
00604
00605 if (depth==2 || abs(ieta)==29 )
00606 hfPET_->HFSetFlagFromPET(*i,*rec,myqual,mySeverity);
00607 }
00608
00609
00610 for (HFRecHitCollection::iterator i = rec->begin();i!=rec->end();++i)
00611 {
00612 int depth=i->id().depth();
00613 int ieta=i->id().ieta();
00614
00615 if (depth==2 || abs(ieta)==29 )
00616 hfS8S1_->HFSetFlagFromS9S1(*i,*rec,myqual,mySeverity);
00617 }
00618
00619
00620 for (HFRecHitCollection::iterator i = rec->begin();i!=rec->end();++i)
00621 {
00622 int depth=i->id().depth();
00623 int ieta=i->id().ieta();
00624
00625 if (depth==1 && abs(ieta)!=29 )
00626 hfS9S1_->HFSetFlagFromS9S1(*i,*rec,myqual, mySeverity);
00627 }
00628 }
00629
00630
00631 e.put(rec);
00632 } else if (subdet_==HcalOther && subdetOther_==HcalCalibration) {
00633 edm::Handle<HcalCalibDigiCollection> digi;
00634 e.getByLabel(inputLabel_,digi);
00635
00636
00637 std::auto_ptr<HcalCalibRecHitCollection> rec(new HcalCalibRecHitCollection);
00638 rec->reserve(digi->size());
00639
00640 int first = firstSample_;
00641 int toadd = samplesToAdd_;
00642
00643 HcalCalibDigiCollection::const_iterator i;
00644 for (i=digi->begin(); i!=digi->end(); i++) {
00645 HcalCalibDetId cell = i->id();
00646
00647 DetId detcell=(DetId)cell;
00648
00649 const HcalChannelStatus* mydigistatus=myqual->getValues(detcell.rawId());
00650 if (mySeverity->dropChannel(mydigistatus->getValue() ) ) continue;
00651 if (dropZSmarkedPassed_)
00652 if (i->zsMarkAndPass()) continue;
00653
00654 const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00655 const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00656 HcalCoderDb coder (*channelCoder, *shape);
00657
00658
00659 if(tsFromDB_) {
00660 const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00661 first = param_ts->firstSample();
00662 toadd = param_ts->samplesToAdd();
00663 }
00664 rec->push_back(reco_.reconstruct(*i,first,toadd,coder,calibrations));
00665
00666
00667
00668
00669
00670
00671
00672
00673
00674
00675
00676
00677
00678
00679 }
00680
00681 e.put(rec);
00682 }
00683 }
00684 delete myqual;
00685 }