|
|
(14 intermediate revisions by the same user not shown) |
Line 1: |
Line 1: |
− | == Histogramming ω Invariant Mass in b1pi ==
| |
− | * Assumes the DReaction has been setup as shown in: [[Analysis_DReaction | DReaction]]
| |
| | | |
− | === Minimal/Simplest ===
| |
− | * Minimal: No check against double-counting, using kinematic fit data only.
| |
− | ** Not best practice: intended only as a simple illustration.
| |
− | * In your DEventProcessor class definition (header):
| |
− | <syntaxhighlight>
| |
− | TH1D* dHist_omegaMass;
| |
− | </syntaxhighlight>
| |
− |
| |
− | * In your DEventProcessor init() method:
| |
− | <syntaxhighlight>
| |
− | dHist_omegaMass = new TH1D("OmegaInvariantMass", ";#pi^{+}#pi^{-}#pi^{0} Invariant Mass (GeV/c^{2})", 600, 0.2, 1.4);
| |
− | </syntaxhighlight>
| |
− |
| |
− | * In your DEventProcessor evnt() method:
| |
− | <syntaxhighlight>
| |
− | vector<const DParticleCombo*> locParticleCombos;
| |
− | locEventLoop->Get(locParticleCombos);
| |
− |
| |
− | for(size_t loc_i = 0; loc_i < locParticleCombos.size(); ++loc_i)
| |
− | {
| |
− | if(locParticleCombos[loc_i]->Get_Reaction()->Get_ReactionName() != "b1pi")
| |
− | continue; //wrong DReaction
| |
− |
| |
− | // get omega combo step //0 is photoproduction, 1 is X(2000) decay, 2 is b1 decay, etc.
| |
− | const DParticleComboStep* locParticleComboStep = locParticleCombo->Get_ParticleComboStep(3);
| |
− |
| |
− | // calculate omega p4 with kinematic fit data
| |
− | TLorentzVector locP4(0.0, 0.0, 0.0, 0.0);
| |
− | locP4 += locParticleComboStep->Get_FinalParticle(0)->lorentzMomentum(); //pi+
| |
− | locP4 += locParticleComboStep->Get_FinalParticle(1)->lorentzMomentum(); //pi-
| |
− | locP4 += locParticleComboStep->Get_FinalParticle(2)->lorentzMomentum(); //pi0
| |
− |
| |
− | Get_Application()->RootWriteLock(); //ACQUIRE ROOT LOCK!!
| |
− | {
| |
− | dHist_omegaMass->Fill(locP4.M());
| |
− | }
| |
− | Get_Application()->RootUnLock(); //RELEASE ROOT LOCK!!
| |
− | break;
| |
− | }
| |
− | </syntaxhighlight>
| |
− |
| |
− | === Thorough: Using an Analysis Action ===
| |
− | * Checks against double-counting, shows/allows histogramming measured or kinfit data.
| |
− | * Action class definition (put in a new header file in your plugin):
| |
− | <syntaxhighlight>
| |
− | class DHistogramAction_omegaMass : public DAnalysisAction
| |
− | {
| |
− | //assumes b1pi reaction
| |
− | public:
| |
− | DHistogramAction_omegaMass(const DReaction* locReaction, bool locUseKinFitResultsFlag, string locActionUniqueString = "") :
| |
− | DAnalysisAction(locReaction, "Hist_omegaMass", locUseKinFitResultsFlag, locActionUniqueString) {} //"Hist_omegaMass" should be unique to this class
| |
− |
| |
− | private:
| |
− | bool Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo);
| |
− | void Initialize(JEventLoop* locEventLoop);
| |
− |
| |
− | TH1D* dHist_omegaMass;
| |
− | };
| |
− | </syntaxhighlight>
| |
− |
| |
− | * Initialize() method: create histogram
| |
− | <syntaxhighlight>
| |
− | void DHistogramAction_omegaMass::Initialize(JEventLoop* locEventLoop)
| |
− | {
| |
− | //CREATE THE HISTOGRAMS
| |
− | Get_Application()->RootWriteLock(); //ACQUIRE ROOT LOCK!!
| |
− |
| |
− | CreateAndChangeTo_ActionDirectory(); //go to the correct directory
| |
− |
| |
− | string locHistName = "InvariantMass";
| |
− | if(gDirectory->Get(locHistName.c_str()) != NULL) //already created by another thread
| |
− | dHist_omegaMass = static_cast<TH1D*>(gDirectory->Get(locHistName.c_str()));
| |
− | else
| |
− | dHist_omegaMass = new TH1D(locHistName.c_str(), ";#pi^{+}#pi^{-}#pi^{0} Invariant Mass (GeV/c^{2})", 600, 0.2, 1.4);
| |
− |
| |
− | Get_Application()->RootUnLock(); //RELEASE ROOT LOCK!!
| |
− | }
| |
− | </syntaxhighlight>
| |
− |
| |
− | * Perform_Action() method: fill histogram
| |
− | <syntaxhighlight>
| |
− | bool DHistogramAction_omegaMass::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
| |
− | {
| |
− | // get omega combo step //0 is photoproduction, 1 is X(2000) decay, 2 is b1 decay, etc.
| |
− | size_t locOmegaStepIndex = 3;
| |
− | const DParticleComboStep* locParticleComboStep = locParticleCombo->Get_ParticleComboStep(locOmegaStepIndex);
| |
− |
| |
− | // check to make sure not double-counting
| |
− | if(!Get_UseKinFitResultsFlag()) //set in constructor
| |
− | {
| |
− | //measured data (kinematic fit results will be unique for each DParticleCombo)
| |
− | deque<pair<const DParticleCombo*, bool> > locPreviousParticleCombos;
| |
− | Get_PreviousParticleCombos(locPreviousParticleCombos); // DParticleCombos for which this action has been previously executed for this event: compare against these
| |
− | //can either compare manually or by using the DAnalysisUtilities methods
| |
− | if(Get_AnalysisUtilities()->Find_SimilarCombos(locParticleCombo, locOmegaStepIndex, locPreviousParticleCombos))
| |
− | return true; //dupe: already histed! //return true: this isn't a cut, the combo didn't fail
| |
− | }
| |
− |
| |
− | // two ways to calculate invariant mass:
| |
− | TLorentzVector locP4_v1 = Get_AnalysisUtilities()->Calc_FinalStateP4(locParticleCombo, locOmegaStepIndex, Get_UseKinFitResultsFlag());
| |
− | // or manually:
| |
− | TLorentzVector locP4_v2(0.0, 0.0, 0.0, 0.0);
| |
− | if(!Get_UseKinFitResultsFlag())
| |
− | {
| |
− | // measured data
| |
− | locP4_v2 += locParticleComboStep->Get_FinalParticle_Measured(0)->lorentzMomentum(); //pi+
| |
− | locP4_v2 += locParticleComboStep->Get_FinalParticle_Measured(1)->lorentzMomentum(); //pi-
| |
− |
| |
− | // pi0 isn't measured: must get data from gamma, gamma
| |
− | const DParticleComboStep* locPi0Step = locParticleCombo->Get_ParticleComboStep(4); //pi0 decay
| |
− | locP4_v2 += locPi0Step->Get_FinalParticle_Measured(0)->lorentzMomentum(); //gamma
| |
− | locP4_v2 += locPi0Step->Get_FinalParticle_Measured(1)->lorentzMomentum(); //gamma
| |
− | }
| |
− | else
| |
− | {
| |
− | // kinematic fit data
| |
− | locP4_v2 += locParticleComboStep->Get_FinalParticle(0)->lorentzMomentum(); //pi+
| |
− | locP4_v2 += locParticleComboStep->Get_FinalParticle(1)->lorentzMomentum(); //pi-
| |
− | locP4_v2 += locParticleComboStep->Get_FinalParticle(2)->lorentzMomentum(); //pi0
| |
− | }
| |
− |
| |
− | // lock around histogram filling
| |
− | Get_Application()->RootWriteLock(); //ACQUIRE ROOT LOCK!!
| |
− | {
| |
− | dHist_omegaMass->Fill(locP4_v1.M());
| |
− | }
| |
− | Get_Application()->RootUnLock(); //RELEASE ROOT LOCK!!
| |
− |
| |
− | return true; //return true: this isn't a cut, the combo didn't fail
| |
− | }
| |
− | </syntaxhighlight>
| |
− |
| |
− | * Add to your DReaction anywhere in the analysis chain:
| |
− | <syntaxhighlight>
| |
− | locReaction->Add_AnalysisAction(new DHistogramAction_omegaMass(locReaction, false)); //false: measured data
| |
− | locReaction->Add_AnalysisAction(new DHistogramAction_omegaMass(locReaction, true, "KinFit")); //true: kinfit data //"KinFit:" a unique, distinguishing string (could have been anything)
| |
− | </syntaxhighlight>
| |
− |
| |
− | == Kinematic Fit Confidence Level ==
| |
− | * In your DEventProcessor evnt() method:
| |
− |
| |
− | <syntaxhighlight>
| |
− | vector<const DParticleCombo*> locParticleCombos;
| |
− | locEventLoop->Get(locParticleCombos);
| |
− |
| |
− | for(size_t loc_i = 0; loc_i < locParticleCombos.size(); ++loc_i)
| |
− | {
| |
− | if(locParticleCombos[loc_i]->Get_Reaction()->Get_ReactionName() != "b1pi")
| |
− | continue; //wrong DReaction
| |
− | const DKinFitResults* locKinFitResults = locParticleCombos[loc_i]->Get_KinFitResults();
| |
− | double locConfidenceLevel = locKinFitResults->Get_ConfidenceLevel();
| |
− | }
| |
− | </syntaxhighlight>
| |
− |
| |
− | == Getting the Thrown Track Topology ==
| |
− |
| |
− | * There are several ways of getting the thrown track information. The most direct way of getting it is to ask JANA for the DMCThrown particles (see below). The parentid member variable matches the myid member variable that the particle decayed from. If it did not decay from any particle (e.g. was directly photoproduced by the generator) then the parentid member will equal 0.
| |
− |
| |
− | <syntaxhighlight>
| |
− | vector<const DMCThrown*> locMCThrowns;
| |
− | locEventLoop->Get(locMCThrowns);
| |
− | </syntaxhighlight>
| |