ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc
(Generate patch)

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents):
Revision 1.21 by khahn, Sat May 12 03:00:15 2012 UTC vs.
Revision 1.23 by khahn, Sun May 20 19:08:59 2012 UTC

# Line 37 | Line 37 | using namespace std;
37   #include "PFCandidateCol.h"
38   #include "PileupInfoCol.h"
39   #include "PileupEnergyDensity.h"
40 + #include "MCParticle.h"
41   #include "TriggerMask.h"
42   #include "TriggerTable.h"
43   #include "Names.h"
# Line 91 | Line 92 | vector<vector<unsigned> > zzcutvec;
92   map<unsigned,float>       evtrhoMap;
93   vector<string> cutstrs;
94   bool passes_HLT_MC;
95 + vector<bool>   PFnoPUflag;;
96  
97   //
98   // prototypes
# Line 101 | Line 103 | double getPUWeight(unsigned npu);
103   void setEffiencyWeights(EventData&, WeightStruct& );
104   void initRunLumiRangeMap();
105   void initEvtRhoMap(map<unsigned,float> &);
106 + unsigned makePFnoPUArray(mithep::Array<PFCandidate> * fPFCandidates,
107 +                         vector<bool> &pfNoPileUpFlag,
108 +                         mithep::Array<Vertex>      * vtxArr );
109   //----------------------------------------------------------------------------
110  
111  
# Line 268 | Line 273 | int main(int argc, char** argv)
273    mithep::Array<mithep::PileupInfo>           *puArr         = new mithep::Array<mithep::PileupInfo>();
274    mithep::Array<mithep::PileupEnergyDensity>  *puDArr        = new mithep::Array<mithep::PileupEnergyDensity>();
275    mithep::Array<mithep::Track>                *trkArr        = new mithep::Array<mithep::Track>();
276 +  mithep::Array<mithep::MCParticle>           *mcArr         = new mithep::Array<mithep::MCParticle>();
277    mithep::TriggerMask                         *trigMask      = new mithep::TriggerMask();
278    mithep::TriggerTable                        *hltTable      = new mithep::TriggerTable();
279    vector<string>                              *hltTableStrings  = new vector<string>();
# Line 280 | Line 286 | int main(int argc, char** argv)
286    TString fPileupInfoName(Names::gkPileupInfoBrn);
287    TString fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn);
288    TString fTrackName(Names::gkTrackBrn);
289 +  TString fMCParticleName(Names::gkMCPartBrn);
290    TString fTriggerMaskName(Names::gkHltBitBrn);
291    TString fTriggerTableName(Names::gkHltTableBrn);
292    
# Line 288 | Line 295 | int main(int argc, char** argv)
295    chain->SetBranchAddress(fMuonName,        &muonArr);
296    chain->SetBranchAddress(fPrimVtxName,     &vtxArr);
297    chain->SetBranchAddress(fPFCandidateName, &pfArr);
298 <  if( ctrl.mc ) chain->SetBranchAddress(fPileupInfoName,  &puArr);
298 >  if( ctrl.mc ) {
299 >    chain->SetBranchAddress(fPileupInfoName,  &puArr);
300 >    chain->SetBranchAddress(fMCParticleName,  &mcArr);
301 >  }
302    chain->SetBranchAddress(fPileupEnergyDensityName, &puDArr);
303    chain->SetBranchAddress(fTrackName, &trkArr);
304    chain->SetBranchAddress(fTriggerMaskName, &trigMask);
# Line 320 | Line 330 | int main(int argc, char** argv)
330    for(UInt_t ientry=0; ientry<imax; ientry++)
331      {
332        chain->GetEntry(ientry);
333 +
334 +      /*
335 +      bool found_tau=false;
336 +      if( ctrl.mc) {  // && ctrl.debug
337 +        cout << "nparticles: " << mcArr->GetEntries() << endl;
338 +        cout << "-----------------------------------" << endl;
339 +        for( int i=0; i<mcArr->GetEntries(); i++ ) {
340 +          const mithep::MCParticle *mc = mcArr->At(i);
341 +          cout << "i: " << i << "\t"
342 +               << "id: " << mc->AbsPdgId()
343 +               << endl;
344 +          if( mc->Is(MCParticle::kTau) ) {
345 +            found_tau=true;
346 +          }
347 +        }
348 +      }
349 +      cout << "-----------------------------------" << endl;
350 +      cout << "found tau: " << found_tau << endl;
351 +      */
352 +
353 +      // events that fail bc ele overlaps w/ reco mu
354 +      /*
355 +      if( info->EvtNum() != 67315 &&
356 +          info->EvtNum() != 288693 ) continue;
357 +      */
358 +
359 +      // events where they apparently don't veto electrons ...
360 +      /*
361 +      if( info->EvtNum() != 141713 &&
362 +          info->EvtNum() != 178019 &&
363 +          info->EvtNum() != 60823 &&
364 +          info->EvtNum() != 192795 ) continue;
365 +      */
366 +
367        if(!(ientry%1000)) cerr << "entry: " << ientry << endl;
368        setPV( ctrl, vtxArr, vtx);
369 +      //if (!(setPV( ctrl, vtxArr, vtx)) ) continue;
370  
371 <      cout << "origrho: " << puDArr->At(0)->Rho()
327 <           << "\torigrholoweta: " << puDArr->At(0)->RhoLowEta()
328 <           << endl;
371 >      float rho = evtrhoMap[info->EvtNum()];
372  
373 +      //
374 +      // pfNoPU
375 +      //
376 +      //mithep::Array<PFCandidate> pfNoPileUp;
377 +      PFnoPUflag.clear();
378 +      int pfnopu_size = makePFnoPUArray( pfArr, PFnoPUflag, vtxArr );
379 +      assert(pfnopu_size == pfArr->GetEntries());
380  
381        //
382        // data/MC
# Line 335 | Line 385 | int main(int argc, char** argv)
385          //
386          // xsec & PU weights
387          //
388 <        //      weights.w = getWeight(xstab,entrymap,chain)/ctrl.totalMC;
388 >        weights.w = getWeight(xstab,entrymap,chain)/ctrl.totalMC;
389 >        cout << "xsec weight: " << weights.w << endl;
390          int npu = -1;
391          for(int i=0;i<puArr->GetEntries();i++) {
392            if(puArr->At(i)->GetBunchCrossing() == 0)
# Line 399 | Line 450 | int main(int argc, char** argv)
450        passingLeptons.clear();
451        EventData ret4l =
452  
453 +
454 +        // reference
455          /*
456          apply_HZZ4L_reference2_selection(ctrl, info,
457                                           vtx,
458                                           pfArr,
459 +                                         //rho,
460                                           puDArr,
461                                           electronArr,
462                                           &electronReferencePreSelection,
# Line 413 | Line 467 | int main(int argc, char** argv)
467                                           &muonIDPFSelection,
468                                           &muonReferenceIsoSelection);
469          */
416 –        apply_HZZ4L_reference2_selection(ctrl, info,
417 –                              vtx,
418 –                              pfArr,
419 –                              puDArr,
420 –                              electronArr,
421 –                              &electronReferencePreSelection,
422 –                              &electronReferenceIDMVASelection,
423 –                              &electronIsoMVASelection,
424 –                              muonArr,
425 –                              &muonReferencePreSelection,
426 –                              &muonIDPFSelection,
427 –                              &muonIsoMVASelection);
470  
471 +        // evolved
472 +        apply_HZZ4L_reference2_selection(ctrl, info,
473 +                                         vtx,
474 +                                         pfArr,
475 +                                         //rho,
476 +                                         puDArr,
477 +                                         electronArr,
478 +                                         &electronReferencePreSelection,
479 +                                         //&electronPreSelection,
480 +                                         &electronReferenceIDMVASelection,
481 +                                         &electronIsoMVASelection,
482 +                                         muonArr,
483 +                                         &muonReferencePreSelection,
484 +                                         //&muonPreSelection,
485 +                                         &muonIDPFSelection,
486 +                                         &muonIsoMVASelection);
487  
488  
489        if( ctrl.debug ) cout << endl;
# Line 441 | Line 499 | int main(int argc, char** argv)
499          if( abs(ret4l.Z2leptons[0].type) == 13  ) Z2type = 13;
500        }  
501        
502 + #ifdef SYNC
503        cout  << "bits: " << ret4l.status.selectionBits  << "\t"
504              << "Z1: " << Z1type << "\t"
505              << "Z2: " << Z2type << endl;
506        cout << endl;
507 <                      
507 > #endif      
508        
509        if( ret4l.status.pass() ) {
510          fillAngles(ret4l,angles);
# Line 548 | Line 607 | bool setPV(ControlFlags ctrl,
607      if( ctrl.debug ) cout << "vertex :: " << i << "\tntrks: " << pv->NTracks() << endl;
608      
609      // Select best PV for corrected d0; if no PV passing cuts, the first PV in the collection will be used
610 <    //  if(!pv->IsValid()) continue;
610 >    if(!pv->IsValid())                                continue;
611      if(pv->NTracksFit()       < fMinNTracksFit)       continue;
612      if(pv->Ndof()                 < fMinNdof)         continue;
613      if(fabs(pv->Z()) > fMaxAbsZ)                      continue;
# Line 569 | Line 628 | bool setPV(ControlFlags ctrl,
628        if( ctrl.debug) cout << "new PV set, pt : " << best_sumpt << endl;
629      }
630    }
631 <  
631 >  // sync
632    if(!bestPV) bestPV = vtxArr->At(0);
633 +  //if(!bestPV) return false;
634    vtx.SetPosition(bestPV->X(),bestPV->Y(),bestPV->Z());
635    vtx.SetErrors(bestPV->XErr(),bestPV->YErr(),bestPV->ZErr());
636    return true;
# Line 665 | Line 725 | void initEvtRhoMap( map<unsigned, float>
725    unsigned evt;
726    float rho;
727  
728 <  FILE * f = fopen("evtrho.dat", "r");
728 >  cout << "initialzing EvtRhoMap ";
729 >  //FILE * f = fopen("evtrho.dat", "r");
730 >  //  FILE * f = fopen("allrho.cali.uniq.dat", "r");
731 >  FILE * f = fopen("allrhoAA.cali.uniq.dat", "r");
732 >  int lcount=0;
733    while( fscanf( f, "%u:%f", &evt, &rho ) ) {
734 +    if( feof(f) ) break;
735      evtrhoMap[evt] = rho;
736 +    lcount++;
737 +    if( !(lcount%1000) ) cout << "."; flush(cout);
738    }
739 <  
739 >  cout << " done" << endl;
740 > }
741 >
742 > //----------------------------------------------------------------------------
743 > unsigned makePFnoPUArray(mithep::Array<PFCandidate> * fPFCandidates,
744 >                     vector<bool> &pfNoPileUpflag,
745 >                     mithep::Array<Vertex>      * vtxArr )
746 > //----------------------------------------------------------------------------
747 > {
748 >  for(UInt_t i = 0; i < fPFCandidates->GetEntries(); i++) {
749 >    const PFCandidate *pf = fPFCandidates->At(i);
750 >    assert(pf);
751 >
752 >    if(pf->PFType() == PFCandidate::eHadron) {
753 >      if(pf->HasTrackerTrk() &&
754 >         vtxArr->At(0)->HasTrack(pf->TrackerTrk()) &&
755 >         vtxArr->At(0)->TrackWeight(pf->TrackerTrk()) > 0) {
756 >
757 >        pfNoPileUpflag.push_back(1);
758 >
759 >      } else {
760 >
761 >        Bool_t vertexFound = kFALSE;
762 >        const Vertex *closestVtx = 0;
763 >        Double_t dzmin = 10000;
764 >        
765 >        // loop over vertices
766 >        for(UInt_t j = 0; j < vtxArr->GetEntries(); j++) {
767 >          const Vertex *vtx = vtxArr->At(j);
768 >          assert(vtx);
769 >          
770 >          if(pf->HasTrackerTrk() &&
771 >             vtx->HasTrack(pf->TrackerTrk()) &&
772 >             vtx->TrackWeight(pf->TrackerTrk()) > 0) {
773 >            vertexFound = kTRUE;
774 >            closestVtx = vtx;
775 >            break;
776 >          }
777 >          Double_t dz = fabs(pf->SourceVertex().Z() - vtx->Z());
778 >          if(dz < dzmin) {
779 >            closestVtx = vtx;
780 >            dzmin = dz;
781 >          }
782 >        }
783 >
784 >        //      if(fCheckClosestZVertex) {
785 >        if(1) {
786 >          // Fallback: if track is not associated with any vertex,
787 >          // associate it with the vertex closest in z
788 >          if(vertexFound || closestVtx != vtxArr->At(0)) {
789 >            //      pfPileUp->Add(pf);
790 >            pfNoPileUpflag.push_back(0);
791 >          } else {
792 >            pfNoPileUpflag.push_back(1);
793 >          }
794 >        } else {
795 >          if(vertexFound && closestVtx != vtxArr->At(0)) {
796 >            //      pfPileUp->Add(pf);
797 >            pfNoPileUpflag.push_back(0);
798 >          } else {
799 >            //      PFCandidate * pfNoPileUp->AddNew(); // Ridiculous but that's how it is
800 >            pfNoPileUpflag.push_back(1);
801 >          }
802 >        }
803 >      } //hadron & trk stuff
804 >    } else { // hadron
805 >      //      PFCandidate * ptr = pfNoPileUp->AddNew();
806 >      pfNoPileUpflag.push_back(1);
807 >    }
808 >  } // Loop over PF candidates
809 >
810 >  return pfNoPileUpflag.size();
811   }
812 +

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines