Fork me on GitHub

Changeset 267 in svn for trunk/Delphes.cpp


Ignore:
Timestamp:
Feb 13, 2009, 8:38:04 PM (15 years ago)
Author:
severine ovyn
Message:

remove double

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/Delphes.cpp

    r264 r267  
    7777int main(int argc, char *argv[])
    7878{
    79 unsigned int nnnn=0, mmmm=0;
    8079
    8180  int appargc = 2;
    … …  
    292291  //Output file : contents of the analysis object data
    293292  ExRootTreeWriter *treeWriter = new ExRootTreeWriter(outputfilename, "Analysis");
     293  ExRootTreeBranch *branchTauJet = treeWriter->NewBranch("TauJet", TRootTauJet::Class());
    294294  ExRootTreeBranch *branchJet = treeWriter->NewBranch("Jet", TRootJet::Class());
    295   ExRootTreeBranch *branchTauJet = treeWriter->NewBranch("TauJet", TRootTauJet::Class());
    296295  ExRootTreeBranch *branchElectron = treeWriter->NewBranch("Electron", TRootElectron::Class());
    297296  ExRootTreeBranch *branchMuon = treeWriter->NewBranch("Muon", TRootMuon::Class());
    … …  
    380379      input_particles.clear();
    381380
    382 
    383 /*
    384 if(0) { // OLD SMEARING ALGORITHM without energy flow
    385       // 2.1 Loop over all particles in event
    386       itGen.Reset();
    387       while( (particle = (TRootGenParticle*) itGen.Next()) ) {
    388           int pid = abs(particle->PID);
    389           float iPhi=0,iEta=0;
    390 
    391        
    392           // 2.1.1********************* preparation for the b-tagging
    393           //// This subarray is needed for the B-jet algorithm
    394           // optimization for speed : put first PID condition, then ETA condition, then either pt or status
    395           if( (pid <= pB || pid == pGLUON) &&// is it a light quark or a gluon, i.e. is it one of these : u,d,c,s,b,g ?
    396               fabs(particle->Eta) < DET->CEN_max_tracker &&
    397               particle->Status != 1 &&
    398               particle->PT > DET->PT_QUARKS_MIN ) {
    399             NFCentralQ.Add(particle);
    400           }
    401          
    402           // 2.1.2 ********************* central detector: keep only visible particles
    403           // keeps only final particles, visible by the central detector, including the fiducial volume
    404           // the ordering of conditions have been optimised for speed : put first the STATUS condition
    405           if( (particle->Status == 1)    &&
    406               (pid != pNU1) && (pid != pNU2) && (pid != pNU3) &&
    407               (fabs(particle->Eta) < DET->CEN_max_calo_fwd)
    408               )
    409             {
    410               genMomentum.SetPxPyPzE(particle->Px, particle->Py, particle->Pz, particle->E);
    411               genMomentumBfield = genMomentum;
    412 
    413               // 2.1.2.1 ********************* central detector: magnetic field
    414               // genMomentumBfield is then changed with respect to the magnetic field
    415               if(DET->FLAG_bfield==1)  TRACP->Propagation(particle,genMomentumBfield);
    416               float eta=fabs(genMomentumBfield.Eta());
    417              
    418  
    419               // 2.1.2.2 ********************* central detector: smearing (calorimeters, muon chambers)
    420               switch(pid) {
    421                
    422               case pE: // all electrons with eta < DET->MAX_CALO_FWD
    423                 DET->SmearElectron(genMomentumBfield);
    424                 if(genMomentumBfield.E()!=0 && eta < DET->CEN_max_tracker && genMomentumBfield.Pt() > DET->PTCUT_elec){
    425                   electron.push_back(ParticleUtil(genMomentumBfield,particle->PID));
    426                 }
    427                 break; // case pE
    428               case pGAMMA: // all photons with eta < DET->MAX_CALO_FWD
    429                 DET->SmearElectron(genMomentumBfield);
    430                 if(genMomentumBfield.E()!=0 && eta < DET->CEN_max_tracker && genMomentumBfield.Pt() > DET->PTCUT_gamma) {
    431                   gamma.push_back(ParticleUtil(genMomentumBfield,particle->PID));
    432                 }
    433                 break; // case pGAMMA
    434               case pMU: // all muons with eta < DET->MAX_MU
    435                 DET->SmearMu(genMomentumBfield);
    436                 if(genMomentumBfield.E()!=0 && eta < DET->CEN_max_mu &&  genMomentumBfield.Pt() > DET->PTCUT_muon){
    437                   muon.push_back(ParticleUtil(genMomentumBfield,particle->PID));
    438                 }
    439                 break; // case pMU
    440               case pLAMBDA: // all lambdas with eta < DET->MAX_CALO_FWD
    441               case pK0S:    // all K0s with eta < DET->MAX_CALO_FWD
    442                 DET->SmearHadron(genMomentumBfield, 0.7);
    443                 break; // case hadron
    444               default:   // all other final particles with eta < DET->MAX_CALO_FWD
    445                 DET->SmearHadron(genMomentumBfield, 1.0);
    446                 break;
    447               } // 2.1.2.2  switch (pid)
    448              
    449  
    450               // 2.1.2.3 ********************* central detector: calotowers
    451               // all final particles but muons and neutrinos   
    452               // for calorimetric towers and missing PT
    453               int charge=Charge(pid);
    454               if(genMomentumBfield.E() !=0 && pid != pMU) {
    455                 // in case the Bfield is not simulated, checks that charged particles have enough pt to reach the calos
    456                 if ( !DET->FLAG_bfield && charge!=0 && genMomentumBfield.Pt() <= DET->TRACK_ptmin ) {}// particules do not reach calos
    457                 else { // particles reach calos
    458                   // applies the calo segmentation and returns iEta & iPhi
    459                   DET->BinEtaPhi(genMomentumBfield.Phi(), genMomentumBfield.Eta(), iPhi, iEta);
    460                   if(iEta != UNDEFINED && iPhi != UNDEFINED) {
    461                       momentumCaloSegmentation.SetPtEtaPhiE(genMomentumBfield.Pt(),iEta,iPhi,genMomentumBfield.E());
    462                       elementCalo = (TRootCalo*) branchCalo->NewEntry();
    463                       elementCalo->Set(momentumCaloSegmentation);
    464                       PhysicsTower Tower(LorentzVector(momentumCaloSegmentation.Px(),momentumCaloSegmentation.Py(),momentumCaloSegmentation.Pz(),momentumCaloSegmentation.E()));
    465                       towers.push_back(Tower);
    466                   } // if iEta != UNDEFINED
    467                 } // else : when particles reach the calos
    468               } // 2.1.2.3 calotowers
    469 
    470 
    471               // 2.1.2.4 ********************* central detector: tracks
    472               // all final charged particles
    473               if(
    474                  (genMomentumBfield.E()!=0) &&
    475                  (fabs(genMomentumBfield.Eta()) < DET->CEN_max_tracker) &&
    476                  (DET->FLAG_bfield || ( !DET->FLAG_bfield && genMomentum.Pt() > DET->TRACK_ptmin )) &&     
    477                         // if bfield not simulated, pt should be high enough to be taken into account
    478                  ((rand()%100) < DET->TRACK_eff)  &&
    479                  (charge!=0)
    480                  )
    481                 {
    482                   elementTracks = (TRootTracks*) branchTracks->NewEntry();
    483                   elementTracks->Set(genMomentum); // fills px,py,pz,pt,e,eta,phi at vertex
    484                   elementTracks->Etaout = genMomentumBfield.Eta();
    485                   elementTracks->Phiout = genMomentumBfield.Phi();
    486                   // TODO!!! apply a smearing on the position of the origin of the track
    487                   // elementTracks->SetPosition(particle->X,particle->Y,particle->Z);
    488                   // uses the output of the bfield computation : Xout=... Yout=... Zout...
    489                   // TODO!!! elementTrakcs->SetPositionOut(Xout,Yout,Zout);
    490                   TrackCentral.push_back(genMomentum); // tracks at vertex!
    491                 } // 2.1.2.4 tracks
    492              
    493             } // 2.1.2 central detector
    494          
    495           // 2.1.3 ********************* forward detectors: zdc
    496           if(DET->FLAG_vfd==1) {
    497              VFD->ZDC(treeWriter,branchZDC,particle);
    498              VFD->RomanPots(treeWriter,branchRP220,branchFP420,particle);
    499           }
    500 } // IF OLD ALGORITHM (= no energy flow )
    501 
    502 
    503 
    504 
    505 else // IF NEW ALGORITHM with energy flow             
    506 */
    507 {
     381      {
    508382      D_CaloTowerList list_of_active_towers;
    509383      D_CaloTowerList list_of_towers_with_photon; // to speed up the code: will only look in interesting towers for gamma candidates
    … …  
    533407            // 2.1a.2.1 Central solenoidal magnetic field
    534408            TRACP->bfield(particle); // fills in particle->EtaCalo et particle->PhiCalo
    535             particle->SetEtaPhi(particle->EtaCalo,particle->PhiCalo); // ???? check this line
    536409
    537410            // 2.1a.2.2 Filling the calorimetric towers -- includes also forward detectors ?
Note: See TracChangeset for help on using the changeset viewer.