#include "CARF/G3Event/interface/G3EventProxy.h" #include "Tracker/TkLayout/interface/CmsTracker.h" #include "Tracker/TkLayout/interface/FullTracker.h" #include "CommonDet/DetLayout/interface/DetLayer.h" #include "CommonDet/BasicDet/interface/SimDet.h" #include "CommonDet/BasicDet/interface/DetUnit.h" #include "CommonDet/BasicDet/interface/DetType.h" #include "Tracker/SiPixelDet/interface/PixelReadout.h" #include "Tracker/SiPixelDet/interface/PixelDet.h" #include "Tracker/SiPixelDet/interface/PixelDetType.h" #include "Tracker/SiStripDet/interface/SiStripDet.h" #include "Tracker/SiStripDet/interface/SiStripDetType.h" #include #include "Utilities/Configuration/interface/Architecture.h" #include "Utilities/Notification/interface/PackageInitializer.h" #include "Utilities/Notification/interface/Observer.h" #include #include "Tracker/SiPixelDet/interface/TrackerMap.h" #include "Tracker/SiPixelDet/interface/TmModule.h" #include "TH1.h" #include "TH2.h" #include "TFile.h" #include "TObjArray.h" #include "TCanvas.h" #include "TString.h" #include "TPostScript.h" using namespace std; #include // forward declare classes we have only pointers to class G3EventProxy; class BuildMap2 : private Observer { public: BuildMap2(); //!< default constructor ~BuildMap2(); //!< default destructor /// this method will do the user analysis void myAnalysis(G3EventProxy * ev); TrackerMap *tkMap; private: /// don't change the name "upDate" - this method is mandatory. void upDate(G3EventProxy * ev) {if (ev!=0) myAnalysis(ev);} int eventsAnalysed; //!< just to count events that have been analysed int PUeventsUsed; //!< just to count events that have been analysed int runsAnalysed; //!< count the runs int lastrun; //!< to remember the last run analysed TObjArray histos; int nhist; }; // ----------------------------------------------------------------------------- // In the constructor we can initialize. It is automatically called when an // instance of the class is created. // ----------------------------------------------------------------------------- BuildMap2::BuildMap2() { cout << "===========================================================" << endl; cout << "=== Start create new BuildMap ===" << endl; // Initialise observer init(); histos(0); tkMap= new TrackerMap(); // Initialise counters eventsAnalysed = 0; PUeventsUsed = 0; runsAnalysed = 0; lastrun = 0; cout << "=== Done create new BuildMap ==" << endl; cout << "===========================================================" << endl; } BuildMap2::~BuildMap2() { tkMap->print(true); cout << "===========================================================" << endl; cout << "=== Start delete BuildMap ===" << endl; cout << " Number of events analysed: " << eventsAnalysed << endl; cout << " Number of runs analysed: " << runsAnalysed << endl; cout << "=== Done delete BuildMap ===" << endl; cout << "===========================================================" << endl; } void BuildMap2::myAnalysis(G3EventProxy * ev) { int simHitsEv = 0; int digiEv = 0; if(eventsAnalysed==0){ tkMap->build(); nhist=0; string s,s1; const CmsTracker::LayerContainer& theLayerContainer = FullTracker::instance()->allLayers(); CmsTracker::layer_p_iter iLay; int iLayer=0; for (iLay = theLayerContainer.begin(); iLay != theLayerContainer.end(); iLay++) {//loop on layers iLayer++; const CmsTracker::DetContainer theDetContainer = (*iLay)->detUnits(); CmsTracker::det_p_iter iDet; int iModule=0; for (iDet = theDetContainer.begin() ; iDet != theDetContainer.end() ; iDet++) { iModule++; int key=iLayer*10000+iModule; TmModule * mod = OrcaModuleMap::omoduleMap[key]; if(mod !=0 && !mod->notInUse()){ ostringstream outs; outs<idModule+mod->ring*1000+mod->layer*100000; outs1 << histname; s1 = outs1.str(); TH1F* h1 = new TH1F(s.c_str(),s1.c_str(),256,0.,256.); histos.Add(h1); mod->histNumber = nhist; nhist++; } } } } if(eventsAnalysed==990){ TFile file("Tracker.root","RECREATE"); histos.Write(); file.Close(); TFile f("Tracker.root"); char name[10]; TString title; string s; int key; Int_t i; for( i=0;i<16540;i++){ if(i%1000==0)cout << i << endl; sprintf(name,"%d",i); TH1F *h = (TH1F*)f.Get(name); title= h->GetTitle(); key = atoi(title.Data()); tkMap->fill(key,float( h->GetMean())); ostringstream outs; outs<< "Mean Adc "<GetMean()<<" RMS ADC "<GetRMS()<<" ;"<GetEntries()<<" digi in this module"; s = outs.str(); cout << s << endl; tkMap->setText(key,s); } tkMap->print(true); f.Close(); } cout << "===========================================================" << endl; cout << "=== Private analysis of " << eventsAnalysed << " event " ; eventsAnalysed++; // some statistics: count events and runs processed CmsTracker::LayerIterator ilayer; const CmsTracker::LayerContainer& bl = FullTracker::instance()->allLayers(); for ( ilayer = bl.begin(); ilayer != bl.end(); ilayer++) { //loop over layers Det::DetUnitContainer theDets = ilayer->detUnits(); Det::DetUnitContainer::iterator idet; for ( idet = theDets.begin(); idet != theDets.end(); idet++) { //loop on detunits TmModule * mod = TmModuleMap::moduleMap[(*idet)]; int digisSize = 0; if(!((mod->layer<32 && mod->layer>28)||(mod->layer<17 &&mod->layer>12))){ SiStripDet* stripDet = (SiStripDet *) *idet; //OK StripReadout & rd = (StripReadout &) (stripDet)->specificReadout(); if( stripDet !=0) { digisSize = rd.ndigis(); } if(digisSize !=0){ StripReadout::Range digis = rd.digis(); for(StripReadout::DigiContainer::const_iterator di=digis.first; di!=digis.second;di++) { // Loop over digis digiEv++; GlobalPoint p = (*idet)->toGlobal(LocalPoint(0,0,0)); float posX =p.x(); float posY =p.y(); float posZ =p.z(); int adc = di->adc(); int strip = di->channel(); //cout << mod->idModule << " " << mod->ring<<" " << mod->layer << " "<histNumber>0) { TH1F * hist = (TH1F*) histos.At(mod->histNumber); hist->Fill(float(adc)); } }//loop over digis }//if digis } else { PixelDet* pixDet = (PixelDet *) *idet; //OK PixelReadout & rd = (PixelReadout &) (pixDet)->specificReadout(); if( pixDet !=0) { digisSize = rd.ndigis(); } // Get the detector module size int numColumns = pixDet->specificTopology().ncolumns(); int numRows = pixDet->specificTopology().nrows(); if(digisSize !=0){ PixelReadout::Range digis = rd.digis(); for(PixelReadout::DigiContainer::const_iterator di=digis.first; di!=digis.second;di++) { // Loop over digis digiEv++; GlobalPoint p = (*idet)->toGlobal(LocalPoint(0,0,0)); float posX =p.x(); float posY =p.y(); float posZ =p.z(); int adc = di->adc(); int strip = di->channel(); int row = di->row(); int column = di->column(); //cout << mod->idModule << " " << mod->ring<<" " << mod->layer << " "<histNumber>0) { TH1F * hist = (TH1F*) histos.At(mod->histNumber); hist->Fill(float(adc)); } }//loop over digis }//if digis }//if isPixel }// end loop on detunits }// end loop on layers cout << " digi " << digiEv << endl; if (ev->simSignal()->id().runNumber() != lastrun) { lastrun = ev->simSignal()->id().runNumber(); runsAnalysed++; } cout << "===========================================================" << endl; }; static PKBuilder mybuilder("BuildMap2Builder");