Fork me on GitHub

Changeset 837fa70 in git for readers


Ignore:
Timestamp:
Oct 9, 2016, 2:57:18 PM (8 years ago)
Author:
Geonmo <ry840901@…>
Branches:
ImprovedOutputFile, Timing, dual_readout, llp, master
Children:
a7c9002
Parents:
386f526
Message:

Mix packedGenParticles and prunedGenParticle.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • readers/DelphesCMSFWLite.cpp

    r386f526 r837fa70  
    9999    handleParticle.getByLabel(event, "genParticles");
    100100  }
    101   else if (!((handlePackedParticle.getBranchNameFor(event, "packedGenParticles")).empty()))
    102   {
     101  else if (!((handlePackedParticle.getBranchNameFor(event, "packedGenParticles")).empty()) && !((handleParticle.getBranchNameFor(event,"prunedGenParticles")).empty()))
     102  {
     103    handleParticle.getByLabel(event, "prunedGenParticles");
    103104    handlePackedParticle.getByLabel(event, "packedGenParticles");
    104105  }
     
    159160
    160161  pdg = TDatabasePDG::Instance();
    161   if( !isMiniAOD) {
    162     for(itParticle = handleParticle->begin(); itParticle != handleParticle->end(); ++itParticle)
    163     {
    164       vectorCandidate.push_back(&*itParticle);
    165     }
    166 
    167     for(itParticle = handleParticle->begin(); itParticle != handleParticle->end(); ++itParticle)
    168     {
    169       const reco::GenParticle &particle = *itParticle;
    170 
    171       pid = particle.pdgId();
    172       status = particle.status();
    173       px = particle.px(); py = particle.py(); pz = particle.pz(); e = particle.energy(); mass = particle.mass();
    174       x = particle.vx(); y = particle.vy(); z = particle.vz();
    175 
    176       candidate = factory->NewCandidate();
    177 
    178       candidate->PID = pid;
    179       pdgCode = TMath::Abs(candidate->PID);
    180 
    181       candidate->Status = status;
    182 
    183       if(particle.mother())
    184       {
    185         itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.mother());
    186         if(itCandidate != vectorCandidate.end()) candidate->M1 = distance(vectorCandidate.begin(), itCandidate);
    187       }
    188 
    189       itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(0));
    190       if(itCandidate != vectorCandidate.end()) candidate->D1 = distance(vectorCandidate.begin(), itCandidate);
    191 
    192       itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(particle.numberOfDaughters() - 1));
    193       if(itCandidate != vectorCandidate.end()) candidate->D2 = distance(vectorCandidate.begin(), itCandidate);
    194 
    195       pdgParticle = pdg->GetParticle(pid);
    196       candidate->Charge = pdgParticle ? Int_t(pdgParticle->Charge()/3.0) : -999;
    197       candidate->Mass = mass;
    198 
    199       candidate->Momentum.SetPxPyPzE(px, py, pz, e);
    200 
    201       candidate->Position.SetXYZT(x*10.0, y*10.0, z*10.0, 0.0);
    202 
    203       allParticleOutputArray->Add(candidate);
    204 
    205       if(!pdgParticle) continue;
    206 
    207       if(status == 1)
    208       {
    209         stableParticleOutputArray->Add(candidate);
    210       }
    211       else if(pdgCode <= 5 || pdgCode == 21 || pdgCode == 15)
    212       {
    213         partonOutputArray->Add(candidate);
    214       }
    215     }
    216   }
    217   else {
    218     for(itPackedParticle = handlePackedParticle->begin(); itPackedParticle != handlePackedParticle->end(); ++itPackedParticle)
    219     {
    220       vectorCandidate.push_back(&*itPackedParticle);
    221     }
    222 
    223     for(itPackedParticle = handlePackedParticle->begin(); itPackedParticle != handlePackedParticle->end(); ++itPackedParticle)
    224     {
    225       const pat::PackedGenParticle &particle = *itPackedParticle;
    226  
    227       pid = particle.pdgId();
    228       status = particle.status();
    229       px = particle.px(); py = particle.py(); pz = particle.pz(); e = particle.energy(); mass = particle.mass();
    230       x = particle.vx(); y = particle.vy(); z = particle.vz();
    231 
    232       candidate = factory->NewCandidate();
    233 
    234       candidate->PID = pid;
    235       pdgCode = TMath::Abs(candidate->PID);
    236 
    237       candidate->Status = status;
    238 
    239       if(particle.mother(0))
    240       {
    241         itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.mother(0));
    242         if(itCandidate != vectorCandidate.end()) candidate->M1 = distance(vectorCandidate.begin(), itCandidate);
    243       }
    244 
    245       itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(0));
    246       if(itCandidate != vectorCandidate.end()) candidate->D1 = distance(vectorCandidate.begin(), itCandidate);
    247 
    248       itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(particle.numberOfDaughters() - 1));
    249       if(itCandidate != vectorCandidate.end()) candidate->D2 = distance(vectorCandidate.begin(), itCandidate);
    250 
    251       pdgParticle = pdg->GetParticle(pid);
    252       candidate->Charge = pdgParticle ? Int_t(pdgParticle->Charge()/3.0) : -999;
    253       candidate->Mass = mass;
    254 
    255       candidate->Momentum.SetPxPyPzE(px, py, pz, e);
    256 
    257       candidate->Position.SetXYZT(x*10.0, y*10.0, z*10.0, 0.0);
    258 
    259       allParticleOutputArray->Add(candidate);
    260 
    261       if(!pdgParticle) continue;
    262 
    263       if(status == 1)
    264       {
    265         stableParticleOutputArray->Add(candidate);
    266       }
    267       else if(pdgCode <= 5 || pdgCode == 21 || pdgCode == 15)
    268       {
    269         partonOutputArray->Add(candidate);
    270       }
     162
     163  for(itParticle = handleParticle->begin(); itParticle != handleParticle->end(); ++itParticle)
     164  {
     165    const reco::GenParticle &particle = *itParticle;
     166    if( !isMiniAOD || particle.status() != 1 ) vectorCandidate.push_back(&*itParticle);
     167  }
     168
     169  for(itParticle = handleParticle->begin(); itParticle != handleParticle->end(); ++itParticle)
     170  {
     171    const reco::GenParticle &particle = *itParticle;
     172
     173    pid = particle.pdgId();
     174    status = particle.status();
     175    if( isMiniAOD && particle.status() == 1 ) continue;
     176    px = particle.px(); py = particle.py(); pz = particle.pz(); e = particle.energy(); mass = particle.mass();
     177    x = particle.vx(); y = particle.vy(); z = particle.vz();
     178
     179    candidate = factory->NewCandidate();
     180
     181    candidate->PID = pid;
     182    pdgCode = TMath::Abs(candidate->PID);
     183
     184    candidate->Status = status;
     185
     186    if(particle.mother())
     187    {
     188      itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.mother());
     189      if(itCandidate != vectorCandidate.end()) candidate->M1 = distance(vectorCandidate.begin(), itCandidate);
     190    }
     191
     192    itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(0));
     193    if(itCandidate != vectorCandidate.end()) candidate->D1 = distance(vectorCandidate.begin(), itCandidate);
     194
     195    itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(particle.numberOfDaughters() - 1));
     196    if(itCandidate != vectorCandidate.end()) candidate->D2 = distance(vectorCandidate.begin(), itCandidate);
     197
     198    pdgParticle = pdg->GetParticle(pid);
     199    candidate->Charge = pdgParticle ? Int_t(pdgParticle->Charge()/3.0) : -999;
     200    candidate->Mass = mass;
     201
     202    candidate->Momentum.SetPxPyPzE(px, py, pz, e);
     203
     204    candidate->Position.SetXYZT(x*10.0, y*10.0, z*10.0, 0.0);
     205
     206    allParticleOutputArray->Add(candidate);
     207
     208    if(!pdgParticle) continue;
     209
     210    if( status == 1)
     211    {
     212      // Prevent duplicated particle.
     213      if ( !isMiniAOD ) stableParticleOutputArray->Add(candidate);
     214    }
     215    else if(pdgCode <= 5 || pdgCode == 21 || pdgCode == 15)
     216    {
     217      partonOutputArray->Add(candidate);
     218    }
     219  }
     220  if ( !isMiniAOD) return ;
     221  // For MiniAOD sample,
     222  // Only status==1 particles are saved to packedGenParticles.
     223  for(itPackedParticle = handlePackedParticle->begin(); itPackedParticle != handlePackedParticle->end(); ++itPackedParticle)
     224  {
     225    vectorCandidate.push_back(&*itPackedParticle);
     226  }
     227
     228  for(itPackedParticle = handlePackedParticle->begin(); itPackedParticle != handlePackedParticle->end(); ++itPackedParticle)
     229  {
     230    const pat::PackedGenParticle &particle = *itPackedParticle;
     231
     232    pid = particle.pdgId();
     233    status = particle.status();
     234    px = particle.px(); py = particle.py(); pz = particle.pz(); e = particle.energy(); mass = particle.mass();
     235    x = particle.vx(); y = particle.vy(); z = particle.vz();
     236
     237    candidate = factory->NewCandidate();
     238
     239    candidate->PID = pid;
     240    pdgCode = TMath::Abs(candidate->PID);
     241
     242    candidate->Status = status;
     243
     244    if(particle.mother(0))
     245    {
     246      itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.mother(0));
     247      if(itCandidate != vectorCandidate.end()) candidate->M1 = distance(vectorCandidate.begin(), itCandidate);
     248    }
     249
     250    itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(0));
     251    if(itCandidate != vectorCandidate.end()) candidate->D1 = distance(vectorCandidate.begin(), itCandidate);
     252
     253    itCandidate = find(vectorCandidate.begin(), vectorCandidate.end(), particle.daughter(particle.numberOfDaughters() - 1));
     254    if(itCandidate != vectorCandidate.end()) candidate->D2 = distance(vectorCandidate.begin(), itCandidate);
     255
     256    pdgParticle = pdg->GetParticle(pid);
     257    candidate->Charge = pdgParticle ? Int_t(pdgParticle->Charge()/3.0) : -999;
     258    candidate->Mass = mass;
     259
     260    candidate->Momentum.SetPxPyPzE(px, py, pz, e);
     261
     262    candidate->Position.SetXYZT(x*10.0, y*10.0, z*10.0, 0.0);
     263
     264    allParticleOutputArray->Add(candidate);
     265
     266    if(!pdgParticle) continue;
     267
     268    if(status == 1)
     269    {
     270      stableParticleOutputArray->Add(candidate);
    271271    }
    272272  }
Note: See TracChangeset for help on using the changeset viewer.