Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-18 02:24:40

0001 // Original Author:  Steven Won
0002 // Original Author:  Alan Campbell for castor
0003 //         Created:  Fri May  2 15:34:43 CEST 2008
0004 // Written to replace the combination of HcalPedestalAnalyzer and HcalPedestalAnalysis
0005 // This code runs 1000x faster and produces all outputs from a single run
0006 // (ADC, fC in .txt plus an .xml file)
0007 //
0008 
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/ServiceRegistry/interface/Service.h"
0016 #include "FWCore/Utilities/interface/InputTag.h"
0017 
0018 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0019 #include "DataFormats/Common/interface/Handle.h"
0020 #include "DataFormats/HcalDetId/interface/HcalGenericDetId.h"
0021 #include "DataFormats/HcalDetId/interface/HcalElectronicsId.h"
0022 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0023 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0024 #include "DataFormats/DetId/interface/DetId.h"
0025 #include "CondFormats/CastorObjects/interface/CastorPedestals.h"
0026 #include "CondFormats/CastorObjects/interface/CastorPedestalWidths.h"
0027 #include "CondFormats/CastorObjects/interface/CastorQIECoder.h"
0028 #include "CondFormats/CastorObjects/interface/CastorQIEData.h"
0029 #include "CondFormats/CastorObjects/interface/CastorQIEShape.h"
0030 #include "CondFormats/CastorObjects/interface/CastorElectronicsMap.h"
0031 #include "CondFormats/CastorObjects/interface/AllObjects.h"
0032 
0033 #include "CalibFormats/CastorObjects/interface/CastorDbRecord.h"
0034 #include "CalibFormats/CastorObjects/interface/CastorDbService.h"
0035 #include "CalibFormats/CastorObjects/interface/CastorCalibrations.h"
0036 #include "CalibFormats/CastorObjects/interface/CastorCalibrationWidths.h"
0037 
0038 #include "CalibCalorimetry/CastorCalib/interface/CastorDbASCIIIO.h"
0039 #include "TBDataFormats/HcalTBObjects/interface/HcalTBTriggerData.h"
0040 
0041 #include "TProfile.h"
0042 #include "TH1.h"
0043 #include "TH2.h"
0044 #include "TCanvas.h"
0045 #include "TStyle.h"
0046 
0047 #include <cmath>
0048 #include <iostream>
0049 #include <map>
0050 #include <memory>
0051 #include <iomanip>
0052 #include <fstream>
0053 #include <vector>
0054 #include <string>
0055 
0056 namespace edm {
0057   class ParameterSet;
0058   class Event;
0059   class EventSetup;
0060 }  // namespace edm
0061 
0062 struct NewPedBunch {
0063   HcalCastorDetId detid;
0064   bool usedflag;
0065   float cap[4];
0066   float capfc[4];
0067   float sig[4][4];
0068   float sigfc[4][4];
0069   float prod[4][4];
0070   float prodfc[4][4];
0071   int num[4][4];
0072 };
0073 
0074 class CastorPedestalsAnalysis : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0075 public:
0076   //Constructor
0077   CastorPedestalsAnalysis(const edm::ParameterSet& ps);
0078   //Destructor
0079   ~CastorPedestalsAnalysis() override;
0080   //Analysis
0081   void analyze(const edm::Event& event, const edm::EventSetup& eventSetup) override;
0082 
0083 private:
0084   //Container for data, 1 per channel
0085   std::vector<NewPedBunch> Bunches;
0086   //Flag for saving histos
0087   bool hiSaveFlag;
0088   bool dumpXML;
0089   bool verboseflag;
0090   int runnum;
0091   int firstTS;
0092   int lastTS;
0093   std::string ROOTfilename;
0094   std::string pedsADCfilename;
0095   std::string pedsfCfilename;
0096   std::string widthsADCfilename;
0097   std::string widthsfCfilename;
0098   std::string XMLfilename;
0099   std::string XMLtag;
0100   std::string ZSfilename;
0101 
0102   edm::ESGetToken<CastorDbService, CastorDbRecord> tok_cond_;
0103   edm::ESGetToken<CastorElectronicsMap, CastorElectronicsMapRcd> tok_map_;
0104 
0105   TH1F* CASTORMeans;
0106   TH1F* CASTORWidths;
0107 
0108   // TH2F *dephist[4];
0109   TH2F* dephist;
0110 
0111   bool firsttime;
0112 
0113   edm::InputTag castorDigiCollectionTag;
0114 };
0115 
0116 CastorPedestalsAnalysis::CastorPedestalsAnalysis(const edm::ParameterSet& ps)
0117     : castorDigiCollectionTag(ps.getParameter<edm::InputTag>("castorDigiCollectionTag")) {
0118   usesResource(TFileService::kSharedResource);
0119 
0120   hiSaveFlag = ps.getUntrackedParameter<bool>("hiSaveFlag", false);
0121   dumpXML = ps.getUntrackedParameter<bool>("dumpXML", false);
0122   verboseflag = ps.getUntrackedParameter<bool>("verbose", false);
0123   firstTS = ps.getUntrackedParameter<int>("firstTS", 0);
0124   lastTS = ps.getUntrackedParameter<int>("lastTS", 9);
0125   firsttime = true;
0126 
0127   tok_cond_ = esConsumes<CastorDbService, CastorDbRecord>();
0128   tok_map_ = esConsumes<CastorElectronicsMap, CastorElectronicsMapRcd>();
0129 }
0130 
0131 CastorPedestalsAnalysis::~CastorPedestalsAnalysis() {
0132   CastorPedestals* rawPedsItem = new CastorPedestals(true);
0133   CastorPedestalWidths* rawWidthsItem = new CastorPedestalWidths(true);
0134   CastorPedestals* rawPedsItemfc = new CastorPedestals(false);
0135   CastorPedestalWidths* rawWidthsItemfc = new CastorPedestalWidths(false);
0136 
0137   //Calculate pedestal constants
0138   std::cout << "Calculating Pedestal constants...\n";
0139   std::vector<NewPedBunch>::iterator bunch_it;
0140   for (bunch_it = Bunches.begin(); bunch_it != Bunches.end(); ++bunch_it) {
0141     if (bunch_it->usedflag) {
0142       if (verboseflag)
0143         std::cout << "Analyzing channel sector= " << bunch_it->detid.sector()
0144                   << " module = " << bunch_it->detid.module() << std::endl;
0145       //pedestal constant is the mean
0146       bunch_it->cap[0] /= bunch_it->num[0][0];
0147       bunch_it->cap[1] /= bunch_it->num[1][1];
0148       bunch_it->cap[2] /= bunch_it->num[2][2];
0149       bunch_it->cap[3] /= bunch_it->num[3][3];
0150       bunch_it->capfc[0] /= bunch_it->num[0][0];
0151       bunch_it->capfc[1] /= bunch_it->num[1][1];
0152       bunch_it->capfc[2] /= bunch_it->num[2][2];
0153       bunch_it->capfc[3] /= bunch_it->num[3][3];
0154       //widths are the covariance matrix--assumed symmetric
0155       bunch_it->sig[0][0] = (bunch_it->prod[0][0] / bunch_it->num[0][0]) - (bunch_it->cap[0] * bunch_it->cap[0]);
0156       bunch_it->sig[0][1] = (bunch_it->prod[0][1] / bunch_it->num[0][1]) - (bunch_it->cap[0] * bunch_it->cap[1]);
0157       bunch_it->sig[0][2] = (bunch_it->prod[0][2] / bunch_it->num[0][2]) - (bunch_it->cap[0] * bunch_it->cap[2]);
0158       bunch_it->sig[0][3] = (bunch_it->prod[0][3] / bunch_it->num[0][3]) - (bunch_it->cap[0] * bunch_it->cap[3]);
0159       bunch_it->sig[1][0] = (bunch_it->prod[1][0] / bunch_it->num[1][0]) - (bunch_it->cap[1] * bunch_it->cap[0]);
0160       bunch_it->sig[1][1] = (bunch_it->prod[1][1] / bunch_it->num[1][1]) - (bunch_it->cap[1] * bunch_it->cap[1]);
0161       bunch_it->sig[1][2] = (bunch_it->prod[1][2] / bunch_it->num[1][2]) - (bunch_it->cap[1] * bunch_it->cap[2]);
0162       bunch_it->sig[1][3] = (bunch_it->prod[1][3] / bunch_it->num[1][3]) - (bunch_it->cap[1] * bunch_it->cap[3]);
0163       bunch_it->sig[2][0] = (bunch_it->prod[2][0] / bunch_it->num[2][0]) - (bunch_it->cap[2] * bunch_it->cap[0]);
0164       bunch_it->sig[2][1] = (bunch_it->prod[2][1] / bunch_it->num[2][1]) - (bunch_it->cap[2] * bunch_it->cap[1]);
0165       bunch_it->sig[2][2] = (bunch_it->prod[2][2] / bunch_it->num[2][2]) - (bunch_it->cap[2] * bunch_it->cap[2]);
0166       bunch_it->sig[2][3] = (bunch_it->prod[2][3] / bunch_it->num[2][3]) - (bunch_it->cap[2] * bunch_it->cap[3]);
0167       bunch_it->sig[3][0] = (bunch_it->prod[3][0] / bunch_it->num[3][0]) - (bunch_it->cap[3] * bunch_it->cap[0]);
0168       bunch_it->sig[3][1] = (bunch_it->prod[3][1] / bunch_it->num[3][1]) - (bunch_it->cap[3] * bunch_it->cap[1]);
0169       bunch_it->sig[3][2] = (bunch_it->prod[3][2] / bunch_it->num[3][2]) - (bunch_it->cap[3] * bunch_it->cap[2]);
0170       bunch_it->sig[3][3] = (bunch_it->prod[3][3] / bunch_it->num[3][3]) - (bunch_it->cap[3] * bunch_it->cap[3]);
0171 
0172       bunch_it->sigfc[0][0] =
0173           (bunch_it->prodfc[0][0] / bunch_it->num[0][0]) - (bunch_it->capfc[0] * bunch_it->capfc[0]);
0174       bunch_it->sigfc[0][1] =
0175           (bunch_it->prodfc[0][1] / bunch_it->num[0][1]) - (bunch_it->capfc[0] * bunch_it->capfc[1]);
0176       bunch_it->sigfc[0][2] =
0177           (bunch_it->prodfc[0][2] / bunch_it->num[0][2]) - (bunch_it->capfc[0] * bunch_it->capfc[2]);
0178       bunch_it->sigfc[0][3] =
0179           (bunch_it->prodfc[0][3] / bunch_it->num[0][3]) - (bunch_it->capfc[0] * bunch_it->capfc[3]);
0180       bunch_it->sigfc[1][0] =
0181           (bunch_it->prodfc[1][0] / bunch_it->num[1][0]) - (bunch_it->capfc[1] * bunch_it->capfc[0]);
0182       bunch_it->sigfc[1][1] =
0183           (bunch_it->prodfc[1][1] / bunch_it->num[1][1]) - (bunch_it->capfc[1] * bunch_it->capfc[1]);
0184       bunch_it->sigfc[1][2] =
0185           (bunch_it->prodfc[1][2] / bunch_it->num[1][2]) - (bunch_it->capfc[1] * bunch_it->capfc[2]);
0186       bunch_it->sigfc[1][3] =
0187           (bunch_it->prodfc[1][3] / bunch_it->num[1][3]) - (bunch_it->capfc[1] * bunch_it->capfc[3]);
0188       bunch_it->sigfc[2][0] =
0189           (bunch_it->prodfc[2][0] / bunch_it->num[2][0]) - (bunch_it->capfc[2] * bunch_it->capfc[0]);
0190       bunch_it->sigfc[2][1] =
0191           (bunch_it->prodfc[2][1] / bunch_it->num[2][1]) - (bunch_it->capfc[2] * bunch_it->capfc[1]);
0192       bunch_it->sigfc[2][2] =
0193           (bunch_it->prodfc[2][2] / bunch_it->num[2][2]) - (bunch_it->capfc[2] * bunch_it->capfc[2]);
0194       bunch_it->sigfc[2][3] =
0195           (bunch_it->prodfc[2][3] / bunch_it->num[2][3]) - (bunch_it->capfc[2] * bunch_it->capfc[3]);
0196       bunch_it->sigfc[3][0] =
0197           (bunch_it->prodfc[3][0] / bunch_it->num[3][0]) - (bunch_it->capfc[3] * bunch_it->capfc[0]);
0198       bunch_it->sigfc[3][1] =
0199           (bunch_it->prodfc[3][1] / bunch_it->num[3][1]) - (bunch_it->capfc[3] * bunch_it->capfc[1]);
0200       bunch_it->sigfc[3][2] =
0201           (bunch_it->prodfc[3][2] / bunch_it->num[3][2]) - (bunch_it->capfc[3] * bunch_it->capfc[2]);
0202       bunch_it->sigfc[3][3] =
0203           (bunch_it->prodfc[3][3] / bunch_it->num[3][3]) - (bunch_it->capfc[3] * bunch_it->capfc[3]);
0204 
0205       for (int i = 0; i != 3; i++) {
0206         CASTORMeans->Fill(bunch_it->cap[i]);
0207         CASTORWidths->Fill(bunch_it->sig[i][i]);
0208       }
0209 
0210       //if(bunch_it->detid.subdet() == 1){
0211 
0212       int fillphi = bunch_it->detid.sector();
0213       //if (bunch_it->detid.depth()==4) fillphi++;
0214 
0215       //    dephist[bunch_it->detid.module()-1]->Fill(bunch_it->detid.ieta(),fillphi,
0216       //             (bunch_it->cap[0]+bunch_it->cap[1]+bunch_it->cap[2]+bunch_it->cap[3])/4);
0217       dephist->Fill(bunch_it->detid.module(),
0218                     fillphi,
0219                     (bunch_it->cap[0] + bunch_it->cap[1] + bunch_it->cap[2] + bunch_it->cap[3]) / 4);
0220 
0221       const CastorPedestal item(bunch_it->detid,
0222                                 bunch_it->cap[0],
0223                                 bunch_it->cap[1],
0224                                 bunch_it->cap[2],
0225                                 bunch_it->cap[3],
0226                                 bunch_it->sig[0][0],
0227                                 bunch_it->sig[1][1],
0228                                 bunch_it->sig[2][2],
0229                                 bunch_it->sig[3][3]);
0230       rawPedsItem->addValues(item);
0231       CastorPedestalWidth widthsp(bunch_it->detid);
0232       widthsp.setSigma(0, 0, bunch_it->sig[0][0]);
0233       widthsp.setSigma(0, 1, bunch_it->sig[0][1]);
0234       widthsp.setSigma(0, 2, bunch_it->sig[0][2]);
0235       widthsp.setSigma(0, 3, bunch_it->sig[0][3]);
0236       widthsp.setSigma(1, 0, bunch_it->sig[1][0]);
0237       widthsp.setSigma(1, 1, bunch_it->sig[1][1]);
0238       widthsp.setSigma(1, 2, bunch_it->sig[1][2]);
0239       widthsp.setSigma(1, 3, bunch_it->sig[1][3]);
0240       widthsp.setSigma(2, 0, bunch_it->sig[2][0]);
0241       widthsp.setSigma(2, 1, bunch_it->sig[2][1]);
0242       widthsp.setSigma(2, 2, bunch_it->sig[2][2]);
0243       widthsp.setSigma(2, 3, bunch_it->sig[2][3]);
0244       widthsp.setSigma(3, 0, bunch_it->sig[3][0]);
0245       widthsp.setSigma(3, 1, bunch_it->sig[3][1]);
0246       widthsp.setSigma(3, 2, bunch_it->sig[3][2]);
0247       widthsp.setSigma(3, 3, bunch_it->sig[3][3]);
0248       rawWidthsItem->addValues(widthsp);
0249 
0250       const CastorPedestal itemfc(bunch_it->detid,
0251                                   bunch_it->capfc[0],
0252                                   bunch_it->capfc[1],
0253                                   bunch_it->capfc[2],
0254                                   bunch_it->capfc[3],
0255                                   bunch_it->sigfc[0][0],
0256                                   bunch_it->sigfc[1][1],
0257                                   bunch_it->sigfc[2][2],
0258                                   bunch_it->sigfc[3][3]);
0259       rawPedsItemfc->addValues(itemfc);
0260       CastorPedestalWidth widthspfc(bunch_it->detid);
0261       widthspfc.setSigma(0, 0, bunch_it->sigfc[0][0]);
0262       widthspfc.setSigma(0, 1, bunch_it->sigfc[0][1]);
0263       widthspfc.setSigma(0, 2, bunch_it->sigfc[0][2]);
0264       widthspfc.setSigma(0, 3, bunch_it->sigfc[0][3]);
0265       widthspfc.setSigma(1, 0, bunch_it->sigfc[1][0]);
0266       widthspfc.setSigma(1, 1, bunch_it->sigfc[1][1]);
0267       widthspfc.setSigma(1, 2, bunch_it->sigfc[1][2]);
0268       widthspfc.setSigma(1, 3, bunch_it->sigfc[1][3]);
0269       widthspfc.setSigma(2, 0, bunch_it->sigfc[2][0]);
0270       widthspfc.setSigma(2, 1, bunch_it->sigfc[2][1]);
0271       widthspfc.setSigma(2, 2, bunch_it->sigfc[2][2]);
0272       widthspfc.setSigma(2, 3, bunch_it->sigfc[2][3]);
0273       widthspfc.setSigma(3, 0, bunch_it->sigfc[3][0]);
0274       widthspfc.setSigma(3, 1, bunch_it->sigfc[3][1]);
0275       widthspfc.setSigma(3, 2, bunch_it->sigfc[3][2]);
0276       widthspfc.setSigma(3, 3, bunch_it->sigfc[3][3]);
0277       rawWidthsItemfc->addValues(widthspfc);
0278     }
0279   }
0280 
0281   // dump the resulting list of pedestals into a file
0282   std::ofstream outStream1(pedsADCfilename.c_str());
0283   CastorDbASCIIIO::dumpObject(outStream1, (*rawPedsItem));
0284   std::ofstream outStream2(widthsADCfilename.c_str());
0285   CastorDbASCIIIO::dumpObject(outStream2, (*rawWidthsItem));
0286 
0287   std::ofstream outStream3(pedsfCfilename.c_str());
0288   CastorDbASCIIIO::dumpObject(outStream3, (*rawPedsItemfc));
0289   std::ofstream outStream4(widthsfCfilename.c_str());
0290   CastorDbASCIIIO::dumpObject(outStream4, (*rawWidthsItemfc));
0291 
0292   if (dumpXML) {
0293     std::ofstream outStream5(XMLfilename.c_str());
0294     //  CastorCondXML::dumpObject (outStream5, runnum, runnum, runnum, XMLtag, 1, (*rawPedsItem), (*rawWidthsItem));
0295   }
0296 
0297   /*
0298   if (hiSaveFlag) {
0299     theFile->Write();
0300   } else {
0301     theFile->cd();
0302     theFile->cd("CASTOR");
0303     CASTORMeans->Write();
0304     CASTORWidths->Write();
0305   }
0306   theFile->cd();
0307   dephist->Write();
0308   */
0309   dephist->SetDrawOption("colz");
0310   dephist->GetXaxis()->SetTitle("module");
0311   dephist->GetYaxis()->SetTitle("sector");
0312 
0313   //for (int n=0; n!= 4; n++)
0314   //{
0315   //dephist[n]->Write();
0316   //dephist[n]->SetDrawOption("colz");
0317   //dephist[n]->GetXaxis()->SetTitle("i#eta");
0318   //dephist[n]->GetYaxis()->SetTitle("i#phi");
0319   //}
0320 
0321   std::stringstream tempstringout;
0322   tempstringout << runnum;
0323   std::string name1 = tempstringout.str() + "_pedplots_1d.png";
0324   std::string name2 = tempstringout.str() + "_pedplots_2d.png";
0325 
0326   TStyle* theStyle = new TStyle("style", "null");
0327   theStyle->SetPalette(1, nullptr);
0328   theStyle->SetCanvasDefH(1200);  //Height of canvas
0329   theStyle->SetCanvasDefW(1600);  //Width of canvas
0330 
0331   gStyle = theStyle;
0332   /*
0333     TCanvas * c1 = new TCanvas("c1","graph",1);
0334     c1->Divide(2,2);
0335     c1->cd(1);
0336     CASTORMeans->Draw();
0337     c1->SaveAs(name1.c_str());   
0338 
0339     theStyle->SetOptStat("n");
0340     gStyle = theStyle;
0341 
0342     TCanvas * c2 = new TCanvas("c2","graph",1);
0343  //   c2->Divide(2,2);
0344     c2->cd(1);
0345     dephist->Draw();
0346     dephist->SetDrawOption("colz");
0347     //c2->cd(2);
0348     //dephist[1]->Draw();
0349     //dephist[1]->SetDrawOption("colz");
0350     //c2->cd(3);
0351     //dephist[2]->Draw();
0352     //dephist[2]->SetDrawOption("colz");
0353     //c2->cd(4);
0354     //dephist[3]->Draw();
0355     //dephist[3]->SetDrawOption("colz");
0356     c2->SaveAs(name2.c_str());
0357 */
0358   /*
0359   std::cout << "Writing ROOT file... ";
0360   theFile->Close();
0361   std::cout << "ROOT file closed.\n";
0362   */
0363 }
0364 
0365 // ------------ method called to for each event  ------------
0366 void CastorPedestalsAnalysis::analyze(const edm::Event& e, const edm::EventSetup& iSetup) {
0367   edm::Handle<CastorDigiCollection> castor;
0368   e.getByLabel(castorDigiCollectionTag, castor);
0369 
0370   auto conditions = &iSetup.getData(tok_cond_);
0371   const CastorQIEShape* shape = conditions->getCastorShape();
0372 
0373   if (firsttime) {
0374     runnum = e.id().run();
0375     std::string runnum_string;
0376     std::stringstream tempstringout;
0377     tempstringout << runnum;
0378     runnum_string = tempstringout.str();
0379     ROOTfilename = runnum_string + "-peds_ADC.root";
0380     pedsADCfilename = runnum_string + "-peds_ADC.txt";
0381     pedsfCfilename = runnum_string + "-peds_fC.txt";
0382     widthsADCfilename = runnum_string + "-widths_ADC.txt";
0383     widthsfCfilename = runnum_string + "-widths_fC.txt";
0384     XMLfilename = runnum_string + "-peds_ADC_complete.xml";
0385     XMLtag = "Castor_pedestals_" + runnum_string;
0386 
0387     edm::Service<TFileService> fs;
0388     fs->mkdir("CASTOR");
0389     /*
0390     theFile = new TFile(ROOTfilename.c_str(), "RECREATE");
0391     theFile->cd();
0392     // Create sub-directories
0393     theFile->mkdir("CASTOR");
0394     theFile->cd();
0395     */
0396     CASTORMeans = fs->make<TH1F>("All Ped Means CASTOR", "All Ped Means CASTOR", 100, 0, 9);
0397     CASTORWidths = fs->make<TH1F>("All Ped Widths CASTOR", "All Ped Widths CASTOR", 100, 0, 3);
0398 
0399     dephist = fs->make<TH2F>("Pedestals (ADC)", "All Castor", 14, 0., 14.5, 16, .5, 16.5);
0400     // dephist[0] = fs->make<TH2F>("Pedestals (ADC)","Depth 1",89, -44, 44, 72, .5, 72.5);
0401     // dephist[1] = fs->make<TH2F>("Pedestals (ADC)","Depth 2",89, -44, 44, 72, .5, 72.5);
0402     // dephist[2] = fs->make<TH2F>("Pedestals (ADC)","Depth 3",89, -44, 44, 72, .5, 72.5);
0403     // dephist[3] = fs->make<TH2F>("Pedestals (ADC)","Depth 4",89, -44, 44, 72, .5, 72.5);
0404 
0405     const CastorElectronicsMap* myRefEMap = &iSetup.getData(tok_map_);
0406     std::vector<HcalGenericDetId> listEMap = myRefEMap->allPrecisionId();
0407     for (std::vector<HcalGenericDetId>::const_iterator it = listEMap.begin(); it != listEMap.end(); ++it) {
0408       HcalGenericDetId mygenid(it->rawId());
0409       if (mygenid.isHcalCastorDetId()) {
0410         NewPedBunch a;
0411         HcalCastorDetId chanid(mygenid.rawId());
0412         a.detid = chanid;
0413         a.usedflag = false;
0414         std::string type = "CASTOR";
0415         for (int i = 0; i != 4; i++) {
0416           a.cap[i] = 0;
0417           a.capfc[i] = 0;
0418           for (int j = 0; j != 4; j++) {
0419             a.sig[i][j] = 0;
0420             a.sigfc[i][j] = 0;
0421             a.prod[i][j] = 0;
0422             a.prodfc[i][j] = 0;
0423             a.num[i][j] = 0;
0424           }
0425         }
0426         Bunches.push_back(a);
0427       }
0428     }
0429     firsttime = false;
0430   }
0431 
0432   std::vector<NewPedBunch>::iterator bunch_it;
0433 
0434   for (CastorDigiCollection::const_iterator j = castor->begin(); j != castor->end(); ++j) {
0435     const CastorDataFrame digi = (const CastorDataFrame)(*j);
0436     for (bunch_it = Bunches.begin(); bunch_it != Bunches.end(); ++bunch_it)
0437       if (bunch_it->detid.rawId() == digi.id().rawId())
0438         break;
0439     bunch_it->usedflag = true;
0440     for (int ts = firstTS; ts != lastTS + 1; ts++) {
0441       const CastorQIECoder* coder = conditions->getCastorCoder(digi.id().rawId());
0442       bunch_it->num[digi.sample(ts).capid()][digi.sample(ts).capid()] += 1;
0443       bunch_it->cap[digi.sample(ts).capid()] += digi.sample(ts).adc();
0444       double charge1 = coder->charge(*shape, digi.sample(ts).adc(), digi.sample(ts).capid());
0445       bunch_it->capfc[digi.sample(ts).capid()] += charge1;
0446       bunch_it->prod[digi.sample(ts).capid()][digi.sample(ts).capid()] +=
0447           (digi.sample(ts).adc() * digi.sample(ts).adc());
0448       bunch_it->prodfc[digi.sample(ts).capid()][digi.sample(ts).capid()] += charge1 * charge1;
0449       if ((ts + 1 < digi.size()) && (ts + 1 < lastTS)) {
0450         bunch_it->prod[digi.sample(ts).capid()][digi.sample(ts + 1).capid()] +=
0451             digi.sample(ts).adc() * digi.sample(ts + 1).adc();
0452         double charge2 = coder->charge(*shape, digi.sample(ts + 1).adc(), digi.sample(ts + 1).capid());
0453         bunch_it->prodfc[digi.sample(ts).capid()][digi.sample(ts + 1).capid()] += charge1 * charge2;
0454         bunch_it->num[digi.sample(ts).capid()][digi.sample(ts + 1).capid()] += 1;
0455       }
0456       if ((ts + 2 < digi.size()) && (ts + 2 < lastTS)) {
0457         bunch_it->prod[digi.sample(ts).capid()][digi.sample(ts + 2).capid()] +=
0458             digi.sample(ts).adc() * digi.sample(ts + 2).adc();
0459         double charge2 = coder->charge(*shape, digi.sample(ts + 2).adc(), digi.sample(ts + 2).capid());
0460         bunch_it->prodfc[digi.sample(ts).capid()][digi.sample(ts + 2).capid()] += charge1 * charge2;
0461         bunch_it->num[digi.sample(ts).capid()][digi.sample(ts + 2).capid()] += 1;
0462       }
0463       if ((ts + 3 < digi.size()) && (ts + 3 < lastTS)) {
0464         bunch_it->prod[digi.sample(ts).capid()][digi.sample(ts + 3).capid()] +=
0465             digi.sample(ts).adc() * digi.sample(ts + 3).adc();
0466         double charge2 = coder->charge(*shape, digi.sample(ts + 3).adc(), digi.sample(ts + 3).capid());
0467         bunch_it->prodfc[digi.sample(ts).capid()][digi.sample(ts + 3).capid()] += charge1 * charge2;
0468         bunch_it->num[digi.sample(ts).capid()][digi.sample(ts + 3).capid()] += 1;
0469       }
0470     }
0471   }
0472 
0473   //this is the last brace
0474 }
0475 
0476 //define this as a plug-in
0477 DEFINE_FWK_MODULE(CastorPedestalsAnalysis);