#include <FillTruth.h>
Inheritance diagram for Anp::FillTruth:
Public Member Functions | |
| FillTruth () | |
| virtual | ~FillTruth () |
| bool | Run (Record &record, TObject *ptr=0) |
| void | Config (const Registry ®) |
Private Member Functions | |
| const Truth | Fill (const NtpMCTruth &ntptruth) const |
| void | FillStdHep (Record &record, const TClonesArray &array) const |
| void | Fill (Truth &truth, const vector< const NtpTHEvent * > &vec) |
| void | Fill (Truth &truth, const vector< const NtpTHTrack * > &vec) |
| void | Fill (Truth &truth, const vector< const NtpTHShower * > &vec) |
| bool | IsConsistent (const Record &record) const |
Private Attributes | |
| bool | fStdHep |
| bool | fCheck |
| int | fHepKey |
| unsigned int | fNTruth |
| unsigned int | fNStHep |
| unsigned int | fNRecoEvent |
| unsigned int | fNRecoTrack |
| unsigned int | fNRecoShower |
|
|
Definition at line 35 of file FillTruth.cxx. 00036 :fStdHep(true), 00037 fCheck(true), 00038 fHepKey(-1), 00039 fNTruth(0), 00040 fNStHep(0), 00041 fNRecoEvent(0), 00042 fNRecoTrack(0), 00043 fNRecoShower(0) 00044 { 00045 }
|
|
|
Definition at line 48 of file FillTruth.cxx. References fNRecoEvent, fNRecoShower, fNRecoTrack, fNStHep, and fNTruth. 00049 {
00050 cout << " Number of truths " << fNTruth << endl
00051 << " Number of stdheps " << fNStHep << endl
00052 << " Number of event truths " << fNRecoEvent << endl
00053 << " Number of track truths " << fNRecoTrack << endl
00054 << " Number of shower truths " << fNRecoShower << endl;
00055 }
|
|
|
Reimplemented from Anp::AlgStore. Definition at line 206 of file FillTruth.cxx. References fCheck, fHepKey, fStdHep, Registry::KeyExists(), Anp::Read(), and reg. 00207 {
00208 Anp::Read(reg, "FillTruthCheck", fCheck);
00209 Anp::Read(reg, "FillTruthStdHep", fStdHep);
00210
00211 if(reg.KeyExists("PrintConfig"))
00212 {
00213 cout << "FillTruth::Config" << endl
00214 << " StdHep = " << fStdHep << endl
00215 << " Check = " << fCheck << endl
00216 << " HepKey = " << fHepKey << endl;
00217 }
00218 }
|
|
||||||||||||
|
Definition at line 449 of file FillTruth.cxx. References count, fNRecoShower, Anp::Truth::fReco, Anp::Truth::IsChild(), MSG, and Anp::Truth::TruthIndex(). 00450 {
00451 //
00452 // Fill all "reconstructed shower to truth" records for this Truth object
00453 //
00454
00455 unsigned int count = 0;
00456
00457 for(vector<const NtpTHShower *>::const_iterator it = vec.begin(); it != vec.end(); ++it)
00458 {
00459 const NtpTHShower *th_shower = *it;
00460
00461 if(th_shower -> neumc != truth.TruthIndex() && !truth.IsChild(th_shower -> neustdhep))
00462 {
00463 continue;
00464 }
00465
00466 const TruthReco truth_reco(TruthReco::kShower,
00467 th_shower -> index,
00468 th_shower -> neumc,
00469 th_shower -> neustdhep,
00470 th_shower -> purity,
00471 th_shower -> completeall,
00472 th_shower -> completeslc);
00473
00474 //
00475 // Fill event truth record in this Truth
00476 //
00477 truth.fReco.push_back(truth_reco);
00478
00479 ++count;
00480 }
00481
00482 fNRecoShower += count;
00483
00484 MSG("FillAlg", Msg::kDebug) << "MC truth record with index " << truth.TruthIndex()
00485 << " matches "<< count << " reconstructed showers" << endl;
00486 }
|
|
||||||||||||
|
Definition at line 407 of file FillTruth.cxx. References count, fNRecoTrack, Anp::Truth::fReco, Anp::Truth::IsChild(), MSG, and Anp::Truth::TruthIndex(). 00408 {
00409 //
00410 // Fill all "reconstructed track to truth" records for this Truth object
00411 //
00412
00413 unsigned int count = 0;
00414
00415 for(vector<const NtpTHTrack *>::const_iterator it = vec.begin(); it != vec.end(); ++it)
00416 {
00417 const NtpTHTrack *th_track = *it;
00418
00419 if(th_track -> neumc != truth.TruthIndex() &&
00420 !truth.IsChild(th_track -> neustdhep) &&
00421 !truth.IsChild(th_track -> trkstdhep))
00422 {
00423 continue;
00424 }
00425
00426 const TruthReco truth_reco(TruthReco::kTrack,
00427 th_track -> index,
00428 th_track -> neumc,
00429 th_track -> trkstdhep,
00430 th_track -> purity,
00431 th_track -> completeall,
00432 th_track -> completeslc);
00433
00434 //
00435 // Fill event truth record in this Truth
00436 //
00437 truth.fReco.push_back(truth_reco);
00438
00439 ++count;
00440 }
00441
00442 fNRecoTrack += count;
00443
00444 MSG("FillAlg", Msg::kDebug) << "MC truth record with index " << truth.TruthIndex()
00445 << " matches "<< count << " reconstructed tracks" << endl;
00446 }
|
|
||||||||||||
|
Definition at line 367 of file FillTruth.cxx. References count, fNRecoEvent, Anp::Truth::fReco, Anp::Truth::IsChild(), MSG, and Anp::Truth::TruthIndex(). 00368 {
00369 //
00370 // Fill all "reconstructed event to truth" records for this Truth object
00371 //
00372
00373 unsigned int count = 0;
00374
00375 for(vector<const NtpTHEvent *>::const_iterator it = vec.begin(); it != vec.end(); ++it)
00376 {
00377 const NtpTHEvent *th_event = *it;
00378
00379 if(th_event -> neumc != truth.TruthIndex() && !truth.IsChild(th_event -> neustdhep))
00380 {
00381 continue;
00382 }
00383
00384 const TruthReco truth_reco(TruthReco::kEvent,
00385 th_event -> index,
00386 th_event -> neumc,
00387 th_event -> neustdhep,
00388 th_event -> purity,
00389 th_event -> completeall,
00390 th_event -> completeslc);
00391
00392 //
00393 // Fill event truth record in this Truth
00394 //
00395 truth.fReco.push_back(truth_reco);
00396
00397 ++count;
00398 }
00399
00400 fNRecoEvent += count;
00401
00402 MSG("FillAlg", Msg::kDebug) << "MC truth record with index " << truth.TruthIndex()
00403 << " matches "<< count << " reconstructed events" << endl;
00404 }
|
|
|
||||||||||||
|
Definition at line 309 of file FillTruth.cxx. References Anp::StdHep::child, Anp::Record::GetStdHep(), Anp::StdHep::IdHEP, Anp::StdHep::index, Anp::StdHep::IstHEP, Anp::StdHep::mass, Anp::StdHep::mc_index, MSG, Anp::StdHep::p4, Anp::StdHep::parent, and Anp::StdHep::vtx. Referenced by Run(). 00310 {
00311 //
00312 // Fill all StdHep blocks from Monte-Carlo simulation
00313 //
00314
00315 const int nstdhep = array.GetEntries();
00316
00317 // initiliaze StdHep vector all at once to speed up
00318 // insertion of new elements into vector
00319 std::vector<Anp::StdHep> &svec = record.GetStdHep();
00320 svec.clear();
00321
00322 // initiliaze stdhep elements
00323 svec.insert(svec.begin(), nstdhep, Anp::StdHep());
00324
00325 for(int index = 0; index < nstdhep; ++index)
00326 {
00327 const NtpMCStdHep *stdhep_ptr = dynamic_cast<const NtpMCStdHep *>(array.At(index));
00328 if(!stdhep_ptr)
00329 {
00330 MSG("FillAlg", Msg::kError) << "No stdhep at " << index << " out of " << nstdhep << endl;
00331 continue;
00332 }
00333
00334 if(int(stdhep_ptr -> index) != index)
00335 {
00336 MSG("FillAlg", Msg::kError) << "Wrong StdHep index" << endl;
00337 continue;
00338 }
00339
00340 StdHep &stdhep = svec[index];
00341
00342 stdhep.index = stdhep_ptr -> index;
00343 stdhep.mc_index = stdhep_ptr -> mc;
00344 stdhep.parent[0] = stdhep_ptr -> parent[0];
00345 stdhep.parent[1] = stdhep_ptr -> parent[1];
00346 stdhep.child[0] = stdhep_ptr -> child[0];
00347 stdhep.child[1] = stdhep_ptr -> child[1];
00348 stdhep.IstHEP = stdhep_ptr -> IstHEP;
00349 stdhep.IdHEP = stdhep_ptr -> IdHEP;
00350 stdhep.mass = stdhep_ptr -> mass;
00351
00352 stdhep.p4[0] = stdhep_ptr -> p4[0];
00353 stdhep.p4[1] = stdhep_ptr -> p4[1];
00354 stdhep.p4[2] = stdhep_ptr -> p4[2];
00355 stdhep.p4[3] = stdhep_ptr -> p4[3];
00356
00357 stdhep.vtx[0] = stdhep_ptr -> vtx[0];
00358 stdhep.vtx[1] = stdhep_ptr -> vtx[1];
00359 stdhep.vtx[2] = stdhep_ptr -> vtx[2];
00360 stdhep.vtx[3] = stdhep_ptr -> vtx[3];
00361
00362 ++fNStHep;
00363 }
00364 }
|
|
|
Definition at line 489 of file FillTruth.cxx. References count, Anp::Record::EventBeg(), Anp::Record::EventEnd(), Anp::EventIter, MSG, Lit::Print(), Anp::Record::ShowerBeg(), Anp::Record::ShowerEnd(), Anp::ShowerIter, Anp::Record::StdHepBeg(), Anp::Record::StdHepEnd(), Anp::StdHepIter, Anp::Record::TrackBeg(), Anp::Record::TrackEnd(), Anp::TrackIter, Anp::Record::TruthBeg(), Anp::Record::TruthEnd(), Anp::TruthIter, and Anp::TruthRecoIter. Referenced by Run(). 00490 {
00491 //
00492 // Run a loop over all truth records in a snarl. For each event and track truth record
00493 // check that there exists a corresponding reco event or track. Small fraction of
00494 // tracks are not included with any reco event. Such tracks for now are ignored.
00495 //
00496
00497 for(TruthIter truth = record.TruthBeg(); truth != record.TruthEnd(); ++truth)
00498 {
00499 for(TruthRecoIter reco = truth -> RecoBeg(); reco != truth -> RecoEnd(); ++reco)
00500 {
00501 short found = 0;
00502
00503 if(record.StdHepBeg() != record.StdHepEnd())
00504 {
00505 for(StdHepIter hep = record.StdHepBeg(); hep != record.StdHepEnd(); ++hep)
00506 {
00507 if(reco -> StdHepIndex() == hep -> StdHepIndex())
00508 {
00509 ++found;
00510 }
00511 }
00512
00513 if(found != 1)
00514 {
00515 MSG("FillAlg", Msg::kError) << "Failed to match StdHep to TruthReco" << endl;
00516 reco -> Print();
00517 return false;
00518 }
00519 }
00520
00521 found = std::count(record.EventBeg(), record.EventEnd(), *reco);
00522 if(reco -> IsEvent() && record.EventBeg() != record.EventEnd() && found != 1)
00523 {
00524 MSG("FillAlg", Msg::kError) << "No match for event to TruthReco: found = " << found << endl;
00525 reco -> Print();
00526 return false;
00527 }
00528
00529 found = std::count(record.TrackBeg(), record.TrackEnd(), *reco);
00530 if(reco -> IsTrack() && record.TrackBeg() != record.TrackEnd() && found != 1)
00531 {
00532 MSG("FillAlg", Msg::kError) << "Failed to match track to TruthReco" << endl;
00533 reco -> Print();
00534 return false;
00535 }
00536 }
00537 }
00538
00539 //
00540 // The loop below runs over all reco events and tracks. For each reco event or track
00541 // run a loop over all Truth records. Check that for each reco event or track there
00542 // is only one matching TruthEvent or TruthTrack object.
00543 //
00544
00545 for(EventIter event = record.EventBeg(); event != record.EventEnd(); ++event)
00546 {
00547 short found = 0;
00548 for(TruthIter truth = record.TruthBeg(); truth != record.TruthEnd(); ++truth)
00549 {
00550 //
00551 // check that some Truth record matches event event
00552 //
00553 found += std::count(truth -> RecoBeg(), truth -> RecoEnd(), *event);
00554 }
00555
00556 //
00557 // every event has its truth info filled?
00558 //
00559 if(found != 1)
00560 {
00561 MSG("FillAlg", Msg::kError) << "Failed to find truth info for a reco event" << endl;
00562 return false;
00563 }
00564 }
00565
00566 for(TrackIter track = record.TrackBeg(); track != record.TrackEnd(); ++track)
00567 {
00568 short found = 0;
00569 for(TruthIter truth = record.TruthBeg(); truth != record.TruthEnd(); ++truth)
00570 {
00571 //
00572 // check that some Truth record matches track event
00573 //
00574 found += std::count(truth -> RecoBeg(), truth -> RecoEnd(), *track);
00575 }
00576
00577 //
00578 // every track has its truth info filled?
00579 //
00580 if(found != 1)
00581 {
00582 MSG("FillAlg", Msg::kError) << "Found " << found <<" truth records for reco track" << endl;
00583 return false;
00584 }
00585 }
00586
00587 for(ShowerIter shower = record.ShowerBeg(); shower != record.ShowerEnd(); ++shower)
00588 {
00589 short found = 0;
00590 for(TruthIter truth = record.TruthBeg(); truth != record.TruthEnd(); ++truth)
00591 {
00592 //
00593 // check that some Truth record matches shower event
00594 //
00595 found += std::count(truth -> RecoBeg(), truth -> RecoEnd(), *shower);
00596 }
00597
00598 //
00599 // every shower has its truth info filled?
00600 //
00601 if(found != 1)
00602 {
00603 MSG("FillAlg", Msg::kError) << "Found " << found <<" truth records for reco shower" << endl;
00604 return false;
00605 }
00606 }
00607
00608 return true;
00609 }
|
|
||||||||||||
|
Implements Anp::AlgStore. Definition at line 58 of file FillTruth.cxx. References Anp::Record::Add(), fCheck, fHepKey, FillStdHep(), Anp::Record::FindStdHep(), DataUtil::GetDetector(), Anp::Record::GetNStdHeps(), Anp::Record::GetNTruths(), GetVldContext(), MCEventInfo::had_fs, Anp::Truth::had_fs, infid(), Anp::Truth::init_state, MCEventInfo::initial_state, IsConsistent(), Anp::Truth::isfid, ReweightHelpers::MCEventInfoFilla(), MSG, MCEventInfo::nucleus, Anp::Truth::nucleus, Anp::Record::StdHepEnd(), Anp::StdHepIter, Anp::Record::TrackBegIterator(), Anp::Record::TrackEndIterator(), Anp::TrackIterator, Anp::Truth::TruthIndex(), Anp::Truth::VtxX(), Anp::Truth::VtxY(), and Anp::Truth::VtxZ(). 00059 {
00060 //
00061 // Get and fill Monte-Carlo truth and StdHep information
00062 //
00063
00064 NtpStRecord *ntprec = dynamic_cast<NtpStRecord *>(ptr);
00065 if(!ntprec)
00066 {
00067 const MomNavigator *mom = dynamic_cast<const MomNavigator *> (ptr);
00068 if(mom)
00069 {
00070 ntprec = dynamic_cast<NtpStRecord *>(mom -> GetFragment("NtpStRecord"));
00071 }
00072 else
00073 {
00074 MSG("FillAlg", Msg::kError) << "Failed to find MomNavigator pointer" << endl;
00075 return false;
00076 }
00077 }
00078
00079 if(!ntprec)
00080 {
00081 MSG("FillAlg", Msg::kError) << "Failed to find NtpStRecord pointer" << endl;
00082 return false;
00083 }
00084
00085 //
00086 // Do nothing if running on data events
00087 //
00088 if(ntprec -> GetVldContext() -> GetSimFlag() != SimFlag::kMC)
00089 {
00090 return true;
00091 }
00092
00093 const TClonesArray *stdhep_array = ntprec -> stdhep;
00094 if(!stdhep_array)
00095 {
00096 MSG("FillAlg", Msg::kWarning) << "TClonesArray of NtpMCStdHep doesn't exist" << endl;
00097 return true;
00098 }
00099
00100 //
00101 // If fStdHep option is set than fill StdHep blocks
00102 //
00103 if(fStdHep)
00104 {
00105 FillStdHep(record, *stdhep_array);
00106 }
00107
00108 const TClonesArray *mc_array = ntprec -> mc;
00109 if(!mc_array)
00110 {
00111 MSG("FillAlg", Msg::kError) << "TClonesArray of NtpMCTruth doesn't exist" << endl;
00112 return false;
00113 }
00114
00115 if(!ntprec -> thevt || !ntprec -> thtrk || !ntprec -> thshw)
00116 {
00117 MSG("FillAlg", Msg::kError) << "One or more TClonesArray of NtpTH doesn't exist" << endl;
00118 return false;
00119 }
00120
00121 const vector<const NtpTHEvent *> thevt_vec = Anp::FillRecoTruth<NtpTHEvent> (ntprec -> thevt);
00122 const vector<const NtpTHTrack *> thtrk_vec = Anp::FillRecoTruth<NtpTHTrack> (ntprec -> thtrk);
00123 const vector<const NtpTHShower *> thshw_vec = Anp::FillRecoTruth<NtpTHShower>(ntprec -> thshw);
00124
00125 MSG("FillAlg", Msg::kDebug)
00126 << endl
00127 << " size of NtpMCTruth array = " << mc_array -> GetEntries() << endl
00128 << " size of NtpTHEvent array = " << thevt_vec.size() << endl
00129 << " size of NtpTHTrack array = " << thtrk_vec.size() << endl
00130 << " size of NtpMCStdHep array = " << stdhep_array -> GetEntries() << endl
00131 << " size of NtpTHShower array = " << thshw_vec.size() << endl;
00132
00133 const int nentries = mc_array -> GetEntries();
00134 for(int index = 0; index < nentries; ++index)
00135 {
00136 const NtpMCTruth *truth_ptr = dynamic_cast<const NtpMCTruth *>(mc_array -> At(index));
00137 if(!truth_ptr)
00138 {
00139 MSG("FillAlg", Msg::kWarning) << "No mc truth at " << index << endl;
00140 continue;
00141 }
00142 if(index != truth_ptr -> index)
00143 {
00144 MSG("FillAlg", Msg::kWarning) << "Wrong array index" << endl;
00145 continue;
00146 }
00147
00148 Truth truth = FillTruth::Fill(*truth_ptr);
00149
00150 FillTruth::Fill(truth, thevt_vec);
00151 FillTruth::Fill(truth, thtrk_vec);
00152 FillTruth::Fill(truth, thshw_vec);
00153
00154 //
00155 // Get parameters for neugen reweighting
00156 //
00157 MCEventInfo mcinfo;
00158 ReweightHelpers::MCEventInfoFilla(&mcinfo, ntprec, truth.TruthIndex());
00159
00160 truth.init_state = mcinfo.initial_state;
00161 truth.nucleus = mcinfo.nucleus;
00162 truth.had_fs = mcinfo.had_fs;
00163
00164 //
00165 // Interaction is in true fiducial volume?
00166 //
00167 truth.isfid = infid(ntprec -> GetVldContext() -> GetDetector(), SimFlag::kMC,
00168 truth.VtxX(), truth.VtxY(), truth.VtxZ());
00169
00170 record.Add(truth);
00171 }
00172
00173 MSG("FillAlg", Msg::kDebug)
00174 << endl
00175 << " size of Truth vector = " << record.GetNTruths() << endl
00176 << " size of StdHep vector = " << record.GetNStdHeps() << endl;
00177
00178 if(fCheck && !IsConsistent(record))
00179 {
00180 return false;
00181 }
00182
00183 //
00184 // Add PDG code for best matching hep block to reconstucted tracks
00185 //
00186 if(fHepKey >= 0)
00187 {
00188 for(TrackIterator itrack = record.TrackBegIterator(); itrack != record.TrackEndIterator(); ++itrack)
00189 {
00190 StdHepIter ihep = record.FindStdHep(*itrack);
00191 if(ihep != record.StdHepEnd())
00192 {
00193 itrack -> Add(fHepKey, ihep -> GetIdHEP());
00194 }
00195 else
00196 {
00197 MSG("FillAlg", Msg::kWarning) << "Failed to find StdHep block" << endl;
00198 }
00199 }
00200 }
00201
00202 return true;
00203 }
|
|
|
Definition at line 58 of file FillTruth.h. |
|
|
Definition at line 60 of file FillTruth.h. |
|
|
Definition at line 64 of file FillTruth.h. Referenced by Fill(), and ~FillTruth(). |
|
|
Definition at line 66 of file FillTruth.h. Referenced by Fill(), and ~FillTruth(). |
|
|
Definition at line 65 of file FillTruth.h. Referenced by Fill(), and ~FillTruth(). |
|
|
Definition at line 63 of file FillTruth.h. Referenced by ~FillTruth(). |
|
|
Definition at line 62 of file FillTruth.h. Referenced by ~FillTruth(). |
|
|
Definition at line 57 of file FillTruth.h. Referenced by Config(). |
1.3.9.1