![]() |
![]() |
#include <RecoTBCalo/EcalHVScan/src/EcalHVScanAnalyzer.h>
Public Member Functions | |
virtual void | analyze (const edm::Event &, const edm::EventSetup &) |
virtual void | beginJob (edm::EventSetup const &) |
EcalHVScanAnalyzer (const edm::ParameterSet &) | |
virtual void | endJob () |
~EcalHVScanAnalyzer () | |
Private Member Functions | |
void | initPNTTMap () |
PNfit | maxPNDiode (const EcalPnDiodeDigi &digi) |
std::pair< int, int > | pnFromTT (const EcalTrigTowerDetId &tt) |
Private Attributes | |
fstream | file |
TH1F | h1d_pnd_ [10] |
TH2F | h2d_anfit_bad_ |
TH2F | h2d_anfit_tot_ |
TH1F | h_ampl_ |
TH1F | h_ampl_xtal_ [85][20] |
TH1F | h_jitt_ |
TH1F | h_jitt_xtal_ [85][20] |
TH1F | h_norm_ampl_xtal_ [85][20] |
TH1F | h_pnd_max_ [10] |
TH1F | h_pnd_t0_ [10] |
std::string | hitCollection_ |
std::string | hitProducer_ |
std::string | pndiodeProducer_ |
std::map< int, std::pair< int, int > > | pnFromTTMap_ |
std::string | rootfile_ |
float | tAmpl [85][20] |
float | tAmplPN1 [85][20] |
float | tAmplPN2 [85][20] |
float | tChi2 [85][20] |
int | tiC [85][20] |
int | tiEta [85][20] |
int | tiPhi [85][20] |
int | tiTT [85][20] |
float | tJitter [85][20] |
int | tnXtal |
TTree * | tree_ |
Implementation: <Notes on="" implementation>="">
Definition at line 51 of file EcalHVScanAnalyzer.h.
EcalHVScanAnalyzer::EcalHVScanAnalyzer | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 43 of file EcalHVScanAnalyzer.cc.
References GenMuonPlsPt100GeV_cfg::cout, lat::endl(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), hitCollection_, hitProducer_, initPNTTMap(), pndiodeProducer_, rootfile_, tAmpl, tAmplPN1, tAmplPN2, tChi2, tiC, tiEta, tiPhi, tiTT, tJitter, and tree_.
00045 { 00046 //now do what ever initialization is needed 00047 rootfile_ = iConfig.getUntrackedParameter<std::string>("rootfile","hvscan.root"); 00048 hitCollection_ = iConfig.getParameter<std::string>("hitCollection"); 00049 hitProducer_ = iConfig.getParameter<std::string>("hitProducer"); 00050 pndiodeProducer_ = iConfig.getParameter<std::string>("pndiodeProducer"); 00051 00052 initPNTTMap(); 00053 00054 00055 std::cout << "EcalHVScanAnalyzer: fetching hitCollection: " << hitCollection_.c_str() 00056 << " produced by " << hitProducer_.c_str() << std::endl; 00057 00058 // initialize the tree 00059 tree_ = new TTree("hvscan","HV Scan Analysis"); 00060 tree_->Branch("ampl",tAmpl,"ampl[85][20]/F"); 00061 tree_->Branch("jitter",tJitter,"jitter[85][20]/F"); 00062 tree_->Branch("chi2",tChi2,"chi2[85][20]/F"); 00063 tree_->Branch("PN1",tAmplPN1,"PN1[85][20]/F"); 00064 tree_->Branch("PN2",tAmplPN2,"PN2[85][20]/F"); 00065 tree_->Branch("iC",tiC,"iC[85][20]/I"); 00066 tree_->Branch("iEta",tiEta,"iEta[85][20]/I"); 00067 tree_->Branch("iPhi",tiPhi,"iPhi[85][20]/I"); 00068 tree_->Branch("iTT",tiTT,"iTT[85][20]/I"); 00069 00070 }
EcalHVScanAnalyzer::~EcalHVScanAnalyzer | ( | ) |
Definition at line 74 of file EcalHVScanAnalyzer.cc.
References tree_.
00076 { 00077 // do anything here that needs to be done at desctruction time 00078 // (e.g. close files, deallocate resources etc.) 00079 delete tree_; 00080 }
void EcalHVScanAnalyzer::analyze | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) | [virtual] |
Implements edm::EDAnalyzer.
Definition at line 213 of file EcalHVScanAnalyzer.cc.
References edm::SortedCollection< T, SORT >::begin(), TestMuL1L2Filter_cff::cerr, GenMuonPlsPt100GeV_cfg::cout, edm::SortedCollection< T, SORT >::end(), lat::endl(), exception, first, edm::Event::getByLabel(), h1d_pnd_, h2d_anfit_bad_, h2d_anfit_tot_, h_ampl_, h_ampl_xtal_, h_jitt_, h_jitt_xtal_, h_norm_ampl_xtal_, h_pnd_max_, h_pnd_t0_, hitCollection_, hitProducer_, i, j, PNfit::max, maxPNDiode(), pndiodeProducer_, pnFromTT(), res, edm::second(), edm::SortedCollection< T, SORT >::size(), PNfit::t0, tAmpl, tAmplPN1, tAmplPN2, tChi2, tiC, tiEta, tiPhi, tiTT, tJitter, and tree_.
00213 { 00214 //======================================================================== 00215 00216 using namespace edm; 00217 using namespace cms; 00218 00219 // take a look at PNdiodes 00220 Handle<EcalPnDiodeDigiCollection> h_pndiodes; 00221 try { 00222 iEvent.getByLabel( pndiodeProducer_,h_pndiodes); 00223 } catch ( std::exception& ex ) { 00224 std::cerr << "Error! can't get the EcalPnDiodeDigiCollection object " << std::endl; 00225 } 00226 const EcalPnDiodeDigiCollection* pndiodes = h_pndiodes.product(); 00227 std::cout << "length of EcalPnDiodeDigiCollection: " << pndiodes->size() << std::endl; 00228 00229 // find max of PND signal 00230 std::vector<double> maxAmplPN; 00231 for(EcalPnDiodeDigiCollection::const_iterator ipnd=pndiodes->begin(); ipnd!=pndiodes->end(); ++ipnd) { 00232 //std::cout << "PNDiode Id: " << ipnd->id() << "\tsize: " << ipnd->size() << std::endl; 00233 for(int is=0; is < ipnd->size() ; ++is) { 00234 float x1 = h1d_pnd_[ipnd->id().iPnId()-1].GetBinContent(is+1); 00235 h1d_pnd_[ipnd->id().iPnId()-1].SetBinContent(is+1, x1+ipnd->sample(is).adc()); 00236 } 00237 00238 // 1D fit to PN diode to find the max 00239 PNfit res = maxPNDiode( *ipnd ); 00240 maxAmplPN.push_back( res.max ); 00241 h_pnd_max_[ipnd->id().iPnId()-1].Fill( res.max ); 00242 h_pnd_t0_[ipnd->id().iPnId()-1].Fill( res.t0 ); 00243 00244 //if(ipnd->id().iPnId() == 2 ) 00245 //std::cout << "PNDiode " << ipnd->id() << "\t max amplitude: " << res.max << ", t0: " << res.t0 << std::endl; 00246 } 00247 00248 // fetch the digis and compute signal amplitude 00249 Handle<EcalUncalibratedRecHitCollection> phits; 00250 try { 00251 //std::cout << "EcalHVScanAnalyzer::analyze getting product with label: " << digiProducer_.c_str()<< " prodname: " << digiCollection_.c_str() << endl; 00252 iEvent.getByLabel( hitProducer_, hitCollection_,phits); 00253 //iEvent.getByLabel( hitProducer_, phits); 00254 } catch ( std::exception& ex ) { 00255 std::cerr << "Error! can't get the product " << hitCollection_.c_str() << std::endl; 00256 } 00257 00258 // reset tree variables 00259 for(int i=0;i<85;++i) { 00260 for(int j=0;j<20;++j) { 00261 tAmpl[i][j] = 0.; 00262 tJitter[i][j] = 0.; 00263 tChi2[i][j] = 0.; 00264 tAmplPN1[i][j] = 0.; 00265 tAmplPN2[i][j] = 0.; 00266 tAmplPN2[i][j] = 0.; 00267 tiC[i][j] = 0; 00268 tiEta[i][j] = 0; 00269 tiPhi[i][j] = 0; 00270 tiTT[i][j] = 0; 00271 } 00272 } 00273 00274 // loop over hits 00275 const EcalUncalibratedRecHitCollection* hits = phits.product(); // get a ptr to the product 00276 std::cout << "# of EcalUncalibratedRecHits hits: " << hits->size() << std::endl; 00277 for(EcalUncalibratedRecHitCollection::const_iterator ithit = hits->begin(); ithit != hits->end(); ++ithit) { 00278 00279 00280 EBDetId anid(ithit->id()); 00281 if(ithit->chi2()>0.) { // make sure fit has converged 00282 h_ampl_.Fill(ithit->amplitude()); 00283 h_jitt_.Fill(ithit->jitter()); 00284 00285 // std::cout << "det Id: " << anid << "\t xtal # " <<anid.ic() << std::endl; 00286 h_ampl_xtal_[anid.ieta()-1][anid.iphi()-1].Fill(ithit->amplitude()); 00287 h_jitt_xtal_[anid.ieta()-1][anid.iphi()-1].Fill(ithit->jitter()); 00288 00289 // normalized amplitude - use average of 2 PNdiodes for each group of TTs 00290 double averagePNdiode = 0.5*( maxAmplPN[ pnFromTT(anid.tower()).first ] + 00291 maxAmplPN[ pnFromTT(anid.tower()).second ] ); 00292 00293 double normAmpl = ithit->amplitude() / averagePNdiode; 00294 h_norm_ampl_xtal_[anid.ieta()-1][anid.iphi()-1].Fill(normAmpl); 00295 00296 00297 /* 00298 std::cout << anid.tower() << " iTT: " << anid.tower().iTT() 00299 << " PNdiode: " << pnFromTT(anid.tower()).first << pnFromTT(anid.tower()).second 00300 << std::endl; 00301 */ 00302 } 00303 00304 // compute fraction of bad analytic fits 00305 h2d_anfit_tot_.Fill(anid.ieta(),anid.iphi()); 00306 if(ithit->chi2()<0) { 00307 h2d_anfit_bad_.Fill(anid.ieta(),anid.iphi()); 00308 } 00309 00310 // fill the tree arrays 00311 //tAmpl[tnXtal] = ithit->amplitude(); 00312 tAmpl[anid.ietaAbs()-1][anid.iphi()-1] = ithit->amplitude(); 00313 tJitter[anid.ietaAbs()-1][anid.iphi()-1] = ithit->jitter(); 00314 tChi2[anid.ietaAbs()-1][anid.iphi()-1] = ithit->chi2(); 00315 tAmplPN1[anid.ietaAbs()-1][anid.iphi()-1] = maxAmplPN[ pnFromTT(anid.tower()).first ]; 00316 tAmplPN2[anid.ietaAbs()-1][anid.iphi()-1] = maxAmplPN[ pnFromTT(anid.tower()).second ]; 00317 tiC[anid.ietaAbs()-1][anid.iphi()-1] = anid.ic(); 00318 tiTT[anid.ietaAbs()-1][anid.iphi()-1] = anid.tower().iTT(); 00319 tiEta[anid.ietaAbs()-1][anid.iphi()-1] = anid.ietaAbs(); 00320 tiPhi[anid.ietaAbs()-1][anid.iphi()-1] = anid.iphi(); 00321 00322 00323 // debugging printout 00324 if(ithit->chi2()<0 && false) 00325 std::cout << "analytic fit failed! EcalUncalibratedRecHit id: " 00326 << EBDetId(ithit->id()) << "\t" 00327 << "amplitude: " << ithit->amplitude() << ", jitter: " << ithit->jitter() 00328 << std::endl; 00329 00330 /* 00331 std::cout << "EBDetId: " << anid 00332 << " iCtal: " << anid.ic() 00333 //<< " trig tow: " << ( (anid.tower_ieta()-1)*4 + anid.tower_iphi()) 00334 << "\tTrig tower ID: " << anid.tower() 00335 << "\tiTT: " << anid.tower().iTT() 00336 << "\tTT hindex: " << anid.tower().hashedIndex() 00337 << std::endl; 00338 */ 00339 00340 00341 00342 }//end of loop over hits 00343 // fill tree 00344 tree_->Fill(); 00345 00346 }
void EcalHVScanAnalyzer::beginJob | ( | edm::EventSetup const & | ) | [virtual] |
Reimplemented from edm::EDAnalyzer.
Definition at line 84 of file EcalHVScanAnalyzer.cc.
References h1d_pnd_, h2d_anfit_bad_, h2d_anfit_tot_, h_ampl_, h_ampl_xtal_, h_jitt_, h_jitt_xtal_, h_norm_ampl_xtal_, h_pnd_max_, h_pnd_t0_, i, and j.
00084 { 00085 //======================================================================== 00086 h_ampl_=TH1F("amplitude","Fitted Pulse Shape Amplitude",800,300,2700); 00087 h_jitt_=TH1F("jitter","Fitted Pulse Shape Jitter",500,0.,5.); 00088 00089 // histo for each PN diode 00090 for(int i=0; i<10; ++i) { 00091 char htit[100]; 00092 char pndiodename[100]; 00093 sprintf(pndiodename,"average ADC count vs. sample - pndiode_%2d",i+1); 00094 sprintf(htit,"pndiode_%2d",i+1); 00095 h1d_pnd_[i] = TH1F(htit,pndiodename,50,0.5,50.5); 00096 00097 sprintf(htit,"max_pndiode_%2d",i+1); 00098 sprintf(pndiodename,"max fitted amplitude - pndiode_%2d",i+1); 00099 h_pnd_max_[i] = TH1F(htit,pndiodename,1000,700,1050); 00100 00101 sprintf(htit,"t0_pndiode_%2d",i+1); 00102 sprintf(pndiodename,"fitted t0 - pndiode_%2d",i+1); 00103 h_pnd_t0_[i] = TH1F(htit,pndiodename,100,20.5,50.5); 00104 } 00105 00106 // histo for each xtal 00107 for(int i=0; i<85; ++i) { 00108 for(int j=0; j<20; ++j) { 00109 char htit[100]; 00110 char hdesc[100]; 00111 00112 sprintf(htit,"amplitude_xtal_%d_%d",i+1,j+1); 00113 sprintf(hdesc,"fitted amplitudes xtal eta:%d phi:%d",i+1,j+1); 00114 h_ampl_xtal_[i][j] = TH1F(htit,hdesc,2400,300.5,2700.5); 00115 00116 sprintf(htit,"norm_ampl_xtal_%d_%d",i+1,j+1); 00117 sprintf(hdesc,"normalized amplitudes xtal eta:%d phi:%d",i+1,j+1); 00118 h_norm_ampl_xtal_[i][j] = TH1F(htit,hdesc,200,1.2,2.0); 00119 00120 00121 sprintf(htit,"jitter_xtal_%d_%d",i+1,j+1); 00122 sprintf(hdesc,"fitted jitters xtal eta:%d phi:%d",i+1,j+1); 00123 h_jitt_xtal_[i][j] = TH1F(htit,hdesc,500,0.,5.); 00124 00125 /* 00126 sprintf(htit,"alpha_xtal_%d_%d",i+1,j+1); 00127 sprintf(hdesc,"fitted \\alpha xtal eta:%d phi:%d",i+1,j+1); 00128 h_alpha_xtal_[i][j] = TH1F(htit,hdesc,200,1.,3.); 00129 00130 sprintf(htit,"tp_xtal_%d_%d",i+1,j+1); 00131 sprintf(hdesc,"fitted t_{p} xtal eta:%d phi:%d",i+1,j+1); 00132 h_tp_xtal_[i][j] = TH1F(htit,hdesc,200,1.,3.); 00133 */ 00134 00135 } 00136 } 00137 00138 h2d_anfit_tot_ = TH2F("anfit_tot","analytic fit total",85,0.5,85.5,20,0.5,20.5); 00139 h2d_anfit_bad_ = TH2F("anfit_bad","analytic fit failed",85,0.5,85.5,20,0.5,20.5); 00140 00141 }
Reimplemented from edm::EDAnalyzer.
Definition at line 145 of file EcalHVScanAnalyzer.cc.
References f, f1, h1d_pnd_, h2d_anfit_bad_, h2d_anfit_tot_, h_ampl_, h_ampl_xtal_, h_jitt_, h_jitt_xtal_, h_norm_ampl_xtal_, h_pnd_max_, h_pnd_t0_, i, int, j, rootfile_, and tree_.
00145 { 00146 //======================================================================== 00147 00148 TFile f(rootfile_.c_str(),"RECREATE"); 00149 h_ampl_.Write(); 00150 h_jitt_.Write(); 00151 00152 for(int i=0; i<10; ++i) { 00153 h1d_pnd_[i].Write(); 00154 h_pnd_max_[i].Write(); 00155 h_pnd_t0_[i].Write(); 00156 } 00157 00158 TH1F h_norm_ampl("norm_ampl","normalized amplitudes all xtals",110,1.0,2.1); 00159 for(int i=0; i<85; ++i) { 00160 for(int j=0; j<20; ++j) { 00161 //if(h_ampl_xtal_[i][j].GetEntries()>0) h_ampl_xtal_[i][j].Fit("gaus","QL"); 00162 h_ampl_xtal_[i][j].Write(); 00163 00164 h_jitt_xtal_[i][j].Write(); 00165 00166 if(h_norm_ampl_xtal_[i][j].GetEntries()>10.) { 00167 h_norm_ampl_xtal_[i][j].Fit("gaus","QL"); 00168 TF1* f1 = h_norm_ampl_xtal_[i][j].GetFunction("gaus"); 00169 h_norm_ampl.Fill( f1->GetParameter(1)); // mean of distributions for each xtal 00170 } 00171 h_norm_ampl_xtal_[i][j].Write(); 00172 } 00173 } 00174 h_norm_ampl.Fit("gaus","QL"); 00175 h_norm_ampl.Write(); 00176 00177 h2d_anfit_tot_.SetOption("COLZ"); 00178 h2d_anfit_tot_.SetXTitle("\\eta index"); 00179 h2d_anfit_tot_.SetYTitle("\\phi index"); 00180 h2d_anfit_tot_.Write(); 00181 h2d_anfit_bad_.SetOption("COLZ"); 00182 h2d_anfit_bad_.SetXTitle("\\eta index"); 00183 h2d_anfit_bad_.SetYTitle("\\phi index"); 00184 h2d_anfit_bad_.Write(); 00185 00186 00187 00188 TH2F h2d_anfit_failure_("anfit_failure","fraction of failed analytic fits",85,0.5,85.5,20,0.5,20.5); 00189 for(int i=1; i <= (int)h2d_anfit_tot_.GetNbinsX(); ++i) { 00190 for(int j=1; j <= (int)h2d_anfit_tot_.GetNbinsY(); ++j) { 00191 h2d_anfit_failure_.SetBinContent(i,j,0.); 00192 if( h2d_anfit_tot_.GetBinContent(i,j)>0. ) 00193 h2d_anfit_failure_.SetBinContent(i,j, h2d_anfit_bad_.GetBinContent(i,j)/h2d_anfit_tot_.GetBinContent(i,j)); 00194 } 00195 } 00196 h2d_anfit_failure_.SetOption("COLZ"); 00197 h2d_anfit_failure_.SetXTitle("\\eta index"); 00198 h2d_anfit_failure_.SetYTitle("\\phi index"); 00199 h2d_anfit_failure_.Write(); 00200 00201 // store the tree in the output file 00202 tree_->Write(); 00203 00204 f.Close(); 00205 }
void EcalHVScanAnalyzer::initPNTTMap | ( | ) | [private] |
Definition at line 379 of file EcalHVScanAnalyzer.cc.
References pnFromTTMap_, and std.
Referenced by EcalHVScanAnalyzer().
00379 { 00380 //======================================================================== 00381 00382 using namespace std; 00383 00384 // pairs of PN diodes for groups of trigger towers 00385 pair<int,int> pair05,pair16,pair27,pair38,pair49; 00386 00387 pair05.first=0; 00388 pair05.second=5; 00389 00390 pair16.first=1; 00391 pair16.second=6; 00392 00393 pair27.first=2; 00394 pair27.second=7; 00395 00396 pair38.first=3; 00397 pair38.second=8; 00398 00399 pair49.first=4; 00400 pair49.second=9; 00401 00402 pnFromTTMap_[1] = pair05; 00403 pnFromTTMap_[2] = pair05; 00404 pnFromTTMap_[3] = pair05; 00405 pnFromTTMap_[4] = pair05; 00406 pnFromTTMap_[5] = pair16; 00407 pnFromTTMap_[6] = pair16; 00408 pnFromTTMap_[7] = pair16; 00409 pnFromTTMap_[8] = pair16; 00410 pnFromTTMap_[9] = pair16; 00411 pnFromTTMap_[10] = pair16; 00412 pnFromTTMap_[11] = pair16; 00413 pnFromTTMap_[12] = pair16; 00414 pnFromTTMap_[13] = pair16; 00415 pnFromTTMap_[14] = pair16; 00416 pnFromTTMap_[15] = pair16; 00417 pnFromTTMap_[16] = pair16; 00418 pnFromTTMap_[17] = pair16; 00419 pnFromTTMap_[18] = pair16; 00420 pnFromTTMap_[19] = pair16; 00421 pnFromTTMap_[20] = pair16; 00422 pnFromTTMap_[21] = pair27; 00423 pnFromTTMap_[22] = pair27; 00424 pnFromTTMap_[23] = pair27; 00425 pnFromTTMap_[24] = pair27; 00426 pnFromTTMap_[25] = pair27; 00427 pnFromTTMap_[26] = pair27; 00428 pnFromTTMap_[27] = pair27; 00429 pnFromTTMap_[28] = pair27; 00430 pnFromTTMap_[29] = pair27; 00431 pnFromTTMap_[30] = pair27; 00432 pnFromTTMap_[31] = pair27; 00433 pnFromTTMap_[32] = pair27; 00434 pnFromTTMap_[33] = pair27; 00435 pnFromTTMap_[34] = pair27; 00436 pnFromTTMap_[35] = pair27; 00437 pnFromTTMap_[36] = pair27; 00438 pnFromTTMap_[37] = pair38; 00439 pnFromTTMap_[38] = pair38; 00440 pnFromTTMap_[39] = pair38; 00441 pnFromTTMap_[40] = pair38; 00442 pnFromTTMap_[41] = pair38; 00443 pnFromTTMap_[42] = pair38; 00444 pnFromTTMap_[43] = pair38; 00445 pnFromTTMap_[44] = pair38; 00446 pnFromTTMap_[45] = pair38; 00447 pnFromTTMap_[46] = pair38; 00448 pnFromTTMap_[47] = pair38; 00449 pnFromTTMap_[48] = pair38; 00450 pnFromTTMap_[49] = pair38; 00451 pnFromTTMap_[50] = pair38; 00452 pnFromTTMap_[51] = pair38; 00453 pnFromTTMap_[52] = pair38; 00454 pnFromTTMap_[53] = pair49; 00455 pnFromTTMap_[54] = pair49; 00456 pnFromTTMap_[55] = pair49; 00457 pnFromTTMap_[56] = pair49; 00458 pnFromTTMap_[57] = pair49; 00459 pnFromTTMap_[58] = pair49; 00460 pnFromTTMap_[59] = pair49; 00461 pnFromTTMap_[60] = pair49; 00462 pnFromTTMap_[61] = pair49; 00463 pnFromTTMap_[62] = pair49; 00464 pnFromTTMap_[63] = pair49; 00465 pnFromTTMap_[64] = pair49; 00466 pnFromTTMap_[65] = pair49; 00467 pnFromTTMap_[66] = pair49; 00468 pnFromTTMap_[67] = pair49; 00469 pnFromTTMap_[68] = pair49; 00470 00471 }
PNfit EcalHVScanAnalyzer::maxPNDiode | ( | const EcalPnDiodeDigi & | digi | ) | [private] |
Definition at line 350 of file EcalHVScanAnalyzer.cc.
References EcalFEMSample::adc(), PNfit::max, res, EcalPnDiodeDigi::sample(), and PNfit::t0.
Referenced by analyze().
00350 { 00351 //======================================================================== 00352 const int ns = 50; 00353 double sample[ns],ampl[ns]; 00354 for(int is=0; is < ns ; ++is) { 00355 sample[is] = is; 00356 ampl[is] = pndigi.sample(is).adc(); 00357 } 00358 TGraph gpn(ns,sample,ampl); 00359 TF1 mypol3("mypol3","pol3",25,50); 00360 gpn.Fit("mypol3","QFR","",25,50); 00361 //TF1* pol3 = (TF1*) gpn.GetFunction("mypol3"); 00362 PNfit res; 00363 res.max = mypol3.GetMaximum(); 00364 res.t0 = mypol3.GetMaximumX(); 00365 00366 return res; 00367 }
std::pair< int, int > EcalHVScanAnalyzer::pnFromTT | ( | const EcalTrigTowerDetId & | tt | ) | [private] |
Definition at line 370 of file EcalHVScanAnalyzer.cc.
References EcalTrigTowerDetId::iTT(), and pnFromTTMap_.
Referenced by analyze().
00370 { 00371 //======================================================================== 00372 00373 return pnFromTTMap_[ tt.iTT() ]; 00374 //return pnFromTTMap_[ 1 ]; 00375 00376 }
fstream EcalHVScanAnalyzer::file [private] |
Definition at line 73 of file EcalHVScanAnalyzer.h.
TH1F EcalHVScanAnalyzer::h1d_pnd_[10] [private] |
Definition at line 76 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH2F EcalHVScanAnalyzer::h2d_anfit_bad_ [private] |
Definition at line 86 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH2F EcalHVScanAnalyzer::h2d_anfit_tot_ [private] |
Definition at line 85 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_ampl_ [private] |
Definition at line 74 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_ampl_xtal_[85][20] [private] |
Definition at line 79 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_jitt_ [private] |
Definition at line 75 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_jitt_xtal_[85][20] [private] |
Definition at line 80 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_norm_ampl_xtal_[85][20] [private] |
Definition at line 83 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_pnd_max_[10] [private] |
Definition at line 77 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
TH1F EcalHVScanAnalyzer::h_pnd_t0_[10] [private] |
Definition at line 78 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), beginJob(), and endJob().
std::string EcalHVScanAnalyzer::hitCollection_ [private] |
Definition at line 70 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
std::string EcalHVScanAnalyzer::hitProducer_ [private] |
Definition at line 71 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
std::string EcalHVScanAnalyzer::pndiodeProducer_ [private] |
Definition at line 72 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
std::map<int, std::pair<int,int> > EcalHVScanAnalyzer::pnFromTTMap_ [private] |
std::string EcalHVScanAnalyzer::rootfile_ [private] |
Definition at line 69 of file EcalHVScanAnalyzer.h.
Referenced by EcalHVScanAnalyzer(), and endJob().
float EcalHVScanAnalyzer::tAmpl[85][20] [private] |
Definition at line 92 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
float EcalHVScanAnalyzer::tAmplPN1[85][20] [private] |
Definition at line 92 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
float EcalHVScanAnalyzer::tAmplPN2[85][20] [private] |
Definition at line 92 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
float EcalHVScanAnalyzer::tChi2[85][20] [private] |
Definition at line 92 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
int EcalHVScanAnalyzer::tiC[85][20] [private] |
Definition at line 94 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
int EcalHVScanAnalyzer::tiEta[85][20] [private] |
Definition at line 94 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
int EcalHVScanAnalyzer::tiPhi[85][20] [private] |
Definition at line 94 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
int EcalHVScanAnalyzer::tiTT[85][20] [private] |
Definition at line 94 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
float EcalHVScanAnalyzer::tJitter[85][20] [private] |
Definition at line 92 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), and EcalHVScanAnalyzer().
int EcalHVScanAnalyzer::tnXtal [private] |
Definition at line 90 of file EcalHVScanAnalyzer.h.
TTree* EcalHVScanAnalyzer::tree_ [private] |
Definition at line 89 of file EcalHVScanAnalyzer.h.
Referenced by analyze(), EcalHVScanAnalyzer(), endJob(), and ~EcalHVScanAnalyzer().