Analysis TTreeFormat

From GlueXWiki
Revision as of 01:36, 8 September 2015 by Pmatt (Talk | contribs) (→‎Combo Tracks (If Any))

Jump to: navigation, search

TTree Format - Overview

  • Called the Physics Analysis Root TTree (PART) format.

Data Hierarchy

  • One TTree per DReaction, each stored in a separate ROOT file.
    • e.g., will have different trees for missing/detected recoil proton for the same final state
  • One TTree entry per event.
  • All particle data stored in arrays/TClonesArray's: one array index per particle.
    • Thrown particles
    • Reconstructed neutral/charged hypotheses
    • Combo particle information (including beam photons)
    • Unused beam particles are NOT saved.
  • Event-independent information (e.g. the target, the DReaction decay chain, etc.) is stored in TTree::fUserInfo (a TList*)

Object-Oriented Programming (NOT YET AVAILABLE)

  • Will provide tool to convert tree branch data into the C++ objects for a given event (callable from TSelector::Process())
    • C++ classes (similar to DANA classes): DTreeCombo, DTreeStep, DTreeParticle (compiled into ROOT dictionary & loadable shared library)
    • Can be used to write/execute reusable software that will work for any reaction (e.g. handling double-counting when filling histograms)

TTree Format: Simulated Data

Thrown Non-Particle Data

// EVENT DATA
"RunNumber": UInt_t
"EventNumber": ULong64_t
"MCWeight": Float_t
 
// # PARTICLES //array size of the thrown particle branches
"NumThrown": UInt_t
 
// THROWN REACTION INFO
"NumPIDThrown_FinalState": ULong64_t //the # of thrown final-state particles (+ pi0) of each type (multiplexed in base 10)
                                       //types (in order from 10^0 -> 10^15): g, e+, e-, nu, mu+, mu-, pi0, pi+, pi-, KLong, K+, K-, n, p, p-bar, n-bar
                                       //e.g. particles decaying from final-state particles are NOT included (e.g. photons from pi0, muons from pions, etc.)
                                     //is sum of #-of-PID * 10^ParticleMultiplexPower() (defined in libraries/include/particleType.h)
                                     //ParticleMultiplexPower() returns a different power of 10 for each final-state PID type. 
                                     //A value of 9 should be interpreted as >= 9.  
"PIDThrown_Decaying": ULong64_t //the types of the thrown decaying particles in the event (multiplexed in base 2)
                                //not the quantity of each, just whether or not they were present (1 or 0)
                                //binary power of a PID is given by ParticleMultiplexPower() (defined in libraries/include/particleType.h)
                                //types: most Particle_t's that aren't final state (e.g. lambda, eta, phi, rho0, etc.) see ParticleMultiplexPower()

Thrown Beam Particle

  • All branch names are prefixed with "ThrownBeam__"
//IDENTIFIER
"PID": Int_t //PDG ID value
 
//KINEMATICS: //At the production vertex 
"X4": TLorentzVector
"P4": TLorentzVector

Thrown Products

  • All branch names are prefixed with "Thrown__"
  • NOTE: The only contains particles corresponding to the "FinalState" and "Decaying" tags of DMCThrown.
    • In other words: No resonances, no decay products of final-state particles, and no orphan particles. 
//IDENTIFIERS / MATCHING
"ParentIndex": Int_t["NumThrown"] //the thrown particle array index of the particle this particle decayed from (-1 if none (e.g. photoproduced))
"PID": Int_t["NumThrown"] //PDG ID value
 
//MATCHING //only present if reconstructed data present (i.e. not if thrown-only tree)
"MatchID": Int_t["NumThrown"] //the "ShowerIndex"/"TrackID" of the reconstructed shower/track that it is matched with (-1 for no match)
"MatchFOM": Float_t["NumThrown"] //Neutrals: confidence level //Tracks: #-matched-hits * hit_fraction //(-1 for no match)
 
//KINEMATICS: //Reported at the particle's production vertex 
"X4": TClonesArray(TLorentzVector["NumThrown"])
"P4": TClonesArray(TLorentzVector["NumThrown"])

TTree Format: Combo-Independent Data

Non-Particle Data

// EVENT DATA
"RunNumber": UInt_t
"EventNumber": ULong64_t
 
// # PARTICLES //these are the array sizes for the particle branches
"NumBeam": UInt_t
"NumChargedHypos": UInt_t
"NumNeutralShowers": UInt_t

Beam Particles (If Used in Combo)

  • Only the beam particles that are included in at least one combo are present.
  • All branch names are prefixed with "Beam__"
//ONLY PRESENT IF BEAM USED IN PARTICLE COMBOS
 
//IDENTIFIERS / MATCHING
"PID": Int_t["NumBeam"] //PDG ID value
"IsGenerator": Bool_t["NumBeam"] // kTRUE/kFALSE if matches the generator beam photon (-1 for no match) //only present if simulated data
 
//KINEMATICS: MEASURED //At the target center
"X4_TargetCenter": TClonesArray(TLorentzVector["NumBeam"])
"P4_Measured": TClonesArray(TLorentzVector["NumBeam"])

Charged Track Hypotheses

  • All branch names are prefixed with "ChargedHypo__"
//IDENTIFIERS / MATCHING
"TrackID": Int_t["NumChargedHypos"] //each physical particle has its own # (to keep track of different pid hypotheses for the same particle)
"PID": Int_t["NumChargedHypos"] //PDG ID value
"ThrownIndex": Int_t["NumChargedHypos"] //the array index of the thrown particle it is matched with (-1 for no match) //only present if simulated data
 
//KINEMATICS: MEASURED  //At the production vertex 
"P4_Measured": TClonesArray(TLorentzVector["NumChargedHypos"])
"X4_Measured": TClonesArray(TLorentzVector["NumChargedHypos"]) //t is the measured value in TOF/BCAL/FCAL projected back to Position_Measured
 
//TRACKING INFO:
"NDF_Tracking": UInt_t["NumChargedHypos"]
"ChiSq_Tracking": Float_t["NumChargedHypos"]
"NDF_DCdEdx": UInt_t["NumChargedHypos"]
"ChiSq_DCdEdx": Float_t["NumChargedHypos"]
"dEdx_CDC": Float_t["NumChargedHypos"] //valid for charged only
"dEdx_FDC": Float_t["NumChargedHypos"] //valid for charged only
 
//TIMING INFO
"HitTime": Float_t["NumChargedHypos"] //the system that is hit is in order of preference: BCAL/TOF/FCAL/ST //to determine which, look whether energy was deposited in these systems
"RFDeltaTVar": Float_t["NumChargedHypos"] //Variance of X4_Measured.T() - RFTime, regardless of which RF bunch is chosen. //Can be used to compute timing ChiSq //RF bunch is combo-dependent
 
//HIT ENERGY:
"dEdx_TOF": Float_t["NumChargedHypos"] //valid for charged only
"dEdx_ST": Float_t["NumChargedHypos"] //valid for charged only
"Energy_BCAL": Float_t["NumChargedHypos"]
"Energy_FCAL": Float_t["NumChargedHypos"]
 
//SHOWER MATCHING:
"TrackBCAL_DeltaPhi": Float_t["NumChargedHypos"] //999.0 if not matched
"TrackBCAL_DeltaZ": Float_t["NumChargedHypos"] //999.0 if not matched
"TrackFCAL_DOCA": Float_t["NumChargedHypos"] //999.0 if not matched

Neutral Showers

  • All branch names are prefixed with "NeutralShower__"
  • Discussion on P4 & X4:
    • Note that P4 is not present because it is defined by X4, and X4 is not present because it is defined by the tracks, which are combo-dependent
    • For combo particles, P4 & X4 are listed for each combo
    • If not used in a combo, can be computed using the shower hit information and the vertex & RF-time of your choosing (e.g. combo production-vertex, RF-time)
  • To determine whether is BCAL or FCAL, see which system has non-zero energy
//MATCHING
"ThrownIndex": Int_t["NumNeutralShowers"] //the array index of the thrown particle it is matched with (-1 for no match) //only present if simulated data
 
//SHOWER INFO
"X4_Shower": Float_t["NumNeutralShowers"] //location/time of the reconstructed shower
"Energy_BCAL": Float_t["NumNeutralShowers"] //is 0.0 if shower in FCAL
"Energy_FCAL": Float_t["NumNeutralShowers"] //is 0.0 if shower in BCAL
 
//NEARBY TRACKS
"TrackBCAL_DeltaPhi": Float_t["NumNeutralShowers"] //is delta to nearest track, is 999.0 if no tracks on BCAL
"TrackBCAL_DeltaZ": Float_t["NumNeutralShowers"] //is delta to nearest track, is 999.0 if no tracks on BCAL
"TrackFCAL_DOCA": Float_t["NumNeutralShowers"] //is DOCA to nearest track, is 999.0 if no tracks on FCAL
 
//PHOTON PID INFO
   //Computed using DVertex (best estimate of reaction vertex using all "good" tracks)
   //Can be used to compute timing chisq //is invalid for non-photons, so computed assuming photon PID
"PhotonRFDeltaTVar": Float_t["NumNeutralShowers"] //Variance of DVertexX4.T() - RFTime, regardless of which RF bunch is chosen. //RF bunch is combo-dependent

TTree Format: Combo-Dependent Data

  • All particle combo data is stored in arrays: array entries correspond to different particle combos

Particle-Independent Data

//NUM COMBOS
"NumCombos": Int_t //size of all of the particle-combo-content arrays
 
//RF
"RFTime_Measured": Float_t["NumCombos"]
"RFTime_KinFit": Float_t["NumCombos"] //only if spacetime kinematic fit performed
 
//KINEMATIC FIT
"ChiSq_KinFit": Float_t["NumCombos"] //only if kinematic fit performed
"NDF_KinFit": UInt_t["NumCombos"] //only if kinematic fit performed // = 0 if kinematic fit doesn't converge

Particle Branch-Name Prefixes

Example Reaction (b1pi):

  • γ p →ω, π+, π-, (p)
    • ω → π+, π-, π0
      • π0 → γ γ

Branch Names:

  • Beam: "ComboBeam"
  • Detected: "PiMinus1", "PiPlus1", "PiPlus2", "PiMinus2", "Photon1", "Photon2"
  • Decaying: "DecayingPi0"
  • Missing: "MissingProton"

Combo Beam Particles (If Any)

  • All branch names are prefixed with "ComboBeam__"
    • E.g. "ComboBeam__BeamIndex"
//IDENTIFIER
"BeamIndex": Int_t //each beam particle has its own #
 
//MEASURED INTERACTION VERTEX
"X4_Measured": TClonesArray(TLorentzVector["NumCombos"])
 
//KINEMATICS: KINFIT //At the interaction vertex //only present if kinfit performed
"X4_KinFit": TClonesArray(TLorentzVector["NumCombos"])
"P4_KinFit": TClonesArray(TLorentzVector["NumCombos"])

Combo Tracks (If Any)

  • All branch names are prefixed with the particle name
    • E.g. "Proton__ChargedIndex", "PiMinus1__P4_KinFit"
//IDENTIFIER
"ChargedIndex": Int_t["NumCombos"] //each physical track has its own # (to keep track of different pid hypotheses for the same track)
 
//MEASURED PID INFO
"Beta_Timing_Measured": Float_t["NumCombos"] // = Path_Length/(c*Delta_t)
"ChiSq_Timing_Measured": Float_t["NumCombos"]
"NDF_Timing": UInt_t["NumCombos"] //same for kinfit chisq
 
//KINEMATIC FIT PID INFO //only present if kinfit performed
"Beta_Timing_KinFit": Float_t["NumCombos"] // = Path_Length/(c*Delta_t)
"ChiSq_Timing_KinFit": Float_t["NumCombos"]
 
//KINEMATIC FIT KINEMATICS //only present if kinfit performed
"X4_KinFit": TClonesArray(TLorentzVector["NumCombos"])
"P4_KinFit": TClonesArray(TLorentzVector["NumCombos"])

Combo Neutrals (If Any)

  • All branch names are prefixed with the particle name
    • E.g. "Photon1__ShowerIndex", "Neutron__P4_KinFit"
//IDENTIFIER
"ShowerIndex": Int_t["NumCombos"]
 
//KINEMATICS: MEASURED  //At the production vertex 
"P4_Measured": TClonesArray(TLorentzVector["NumChargedHypos"])
"X4_Measured": TClonesArray(TLorentzVector["NumChargedHypos"]) //t is the measured value in TOF/BCAL/FCAL projected back to Position_Measured
 
//MEASURED PID INFO
"Beta_Timing_Measured": Float_t["NumCombos"] // = Path_Length/(c*Delta_t)
"ChiSq_Timing_Measured": Float_t["NumCombos"] //only present if photon
"NDF_Timing": UInt_t["NumCombos"] //only present if photon //same for kinfit chisq
 
//KINEMATIC FIT PID INFO
"Beta_Timing_KinFit": Float_t["NumCombos"] // = Path_Length/(c*Delta_t)
"ChiSq_Timing_KinFit": Float_t["NumCombos"] //only present if photon
 
//KINEMATIC FIT KINEMATICS //only present if kinfit performed
"X4_KinFit": TClonesArray(TLorentzVector["NumCombos"])
"P4_KinFit": TClonesArray(TLorentzVector["NumCombos"])

Combo Decaying Particles (If Any, If Detached/KinFit)

  • All branch names are prefixed with "Decaying" and the particle name
    • E.g.: "DecayingPi0__X4"
//KINEMATICS: //At the decay vertex 
"X4": TLorentzVector //only present if has a detached vertex //kinematic fit result if kinfit performed, else reconstructed from detected particles
"P4_KinFit": TLorentzVector //only present if kinfit performed

Combo Missing Particles (If Any & If KinFit)

  • All branch names are prefixed with "Missing" and the particle name
    • E.g.: "MissingProton__P4_KinFit"
//KINFIT KINEMATICS: //At its production vertex //only present if kinfit performed
"P4_KinFit": TLorentzVector

TTree Format: DReaction Info

  • Stored in TTree::fUserInfo (a TList*)
  • "ParticleNameList": TList of the names of the reaction particles in the tree, in the order they were specified in the DReaction.
  • "MiscInfoMap": TMap of TObjString -> TObjString
    • "KinFitType" -> DKinFitType (converted to TObjString)
    • "Target__PID" -> int (converted to TObjString): PDG PID of target particle //if a target particle was specified
    • "Target__Mass" -> double (converted to TObjString): Mass of the target particle. //if a target particle was specified
    • "Missing__PID" -> int (converted to TObjString): PDG PID of missing particle //if a missing particle was specified
    • "MissingNAME__Mass" -> double (converted to TObjString): Mass of the 'NAME' missing particle (e.g. 'NAME' = Proton). //if a missing particle was specified
    • "DecayingNAME__Mass" -> double (converted to TObjString): Mass of the 'NAME' decaying particle (e.g. 'NAME' = Pi0). //if decaying particles were present
  • "NameToPIDMap": TMap of "UniqueParticleName" (TObjString) -> int (PDG) (converted to TObjString)
  • "NameToPositionMap": TMap of "UniqueParticleName" (TObjString) -> "StepIndex_ParticleIndex" (stored in TObjString) (ParticleIndex = -1 for initial, -2 for target, 0+ for final state)
  • "PositionToNameMap": TMap of "StepIndex_ParticleIndex" (stored in TObjString) (ParticleIndex = -1 for initial, -2 for target, 0+ for final state) -> "UniqueParticleName" (TObjString)
  • "PositionToPIDMap": TMap of "StepIndex_ParticleIndex" (stored in TObjString) (ParticleIndex = -1 for initial, -2 for target, 0+ for final state) -> int (PDG) (converted to TObjString)
  • "DecayProductMap": TMap of "DecayingParticleName" (TObjString) -> "DecayProductNames" (stored in a TList of TObjString objects). Excludes resonances and intermediate decays (e.g. if Ξ-→π-Λ→π-π-p: will be Ξ-→π-π-p and Λ decay not listed)

Usage

Create TTrees

  • To save data to a TTree for a given DReaction, TTree output must be first be enabled for that reaction. See DReaction Control Variables for details.
    • Note: Only one thrown tree will be created during program execution. If the DEventWriterROOT::Create_ThrownTree() function is called more than once, nothing happens on subsequent calls. 
#include "ANALYSIS/DEventWriterROOT.h"
//In plugin brun():
const DEventWriterROOT* locEventWriterROOT = NULL;
locEventLoop->GetSingle(locEventWriterROOT);
locEventWriterROOT->Create_DataTrees(locEventLoop); //creates TTrees for all output-enabled DReactions
locEventWriterROOT->Create_ThrownTree("tree_b1pi_thrownmc.root"); //optional: create a ttree containing only the thrown data //string is output file name

Save Data to TTree

  • The below only saves the particle combinations (for TTree-output-enabled DReaction's created in the factory specified by the tag) that survived all of the DAnalysisAction cuts. 
//In plugin evnt()
const DEventWriterROOT* locEventWriterROOT = NULL;
locEventLoop->GetSingle(locEventWriterROOT);
locEventWriterROOT->Fill_DataTrees(locEventLoop, "b1pi_hists"); //string is the DReaction factory tag that the DReactions were created in
  • The below allows you to choose which DParticleCombo's (locParticleCombos) of which DReaction's (locReaction) to save.
    • Beware: the locParticleCombos MUST have originated from the locReaction or else this will probably crash (can check DParticleCombo::Get_Reaction()). 
//In plugin evnt()
#include "ANALYSIS/DEventWriterROOT.h"
vector<const DEventWriterROOT*> locEventWriterROOTVector;
locEventLoop->Get(locEventWriterROOTVector); //creates the TTrees for all DReactions upon first call
locEventWriterROOTVector[0]->Fill_Tree(locEventLoop, locReaction, locParticleCombos);
  • The below fills a TTree that only contains the thrown particle data. 
//In plugin evnt()
const DEventWriterROOT* locEventWriterROOT = NULL;
locEventLoop->GetSingle(locEventWriterROOT);
locEventWriterROOT->Fill_ThrownTree(locEventLoop);

Accessing TTree Data

  • TTree:
MyTree->Draw("PiMinus1__P4_Measured->Theta()"); //draws all particle combinations
  • TBrowser (draws all particle combinations):
b1pi Events


GetEntry(entry);
 
//this doesn't take into account double-counting!!
TLorentzVector locPiZeroP4 = *Gamma1__P4_KinFit + *Gamma2__P4_KinFit;
dPiZeroMassHist->Fill(locPiZeroP4.M());
 
TLorentzVector locOmegaP4 = locPiZeroP4 + *PiPlus2__P4_KinFit + *PiMinus2__P4_KinFit;
dOmegaMassHist->Fill(locOmegaP4.M());
 
TLorentzVector locB1PlusP4 = locOmegaP4 + *PiPlus1__P4_KinFit;
dB1PlusMassHist->Fill(locB1PlusP4.M());
 
TLorentzVector locX2000P4 = locB1PlusP4 + *PiMinus1__P4_KinFit;
dX2000MassHist->Fill(locX2000P4.M());

Usage - Advanced

Custom Branches

  • You can create and fill custom branches by inheriting from the DEventWriterROOT class to create your own writer class.
  • Use the trunk/scripts/analysis/MakeEventWriterROOT.pl script to generate the necessary code to do this.
  • Run this perl script with no arguments to get complete usage instructions.

Preventing Double-Counting

  • Since you can have multiple particle combinations per event, you have to be very careful to make sure you aren't double-counting when filling your histograms.
    • For example, if you're histogramming the invariant mass of the π0's decay to γγ in b1pi events using the measured photon data, multiple combinations may use the same showers for the photons, while having different tracks for the other particles. Here's the recommended solution for checking this in your TSelector:
  • Keep track of "ObjectID"'s of the particles for previous particle combinations that pass each cut that you place in your TSelector
    • Put code like this in your header file
//member variables added to your TSelector
vector<vector<int> > sPreviousObjectIDs_Beam; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_PiMinus1; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_PiMinus2; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_PiPlus1; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_PiPlus2; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_Gamma1; //1st dimension is cut-index
vector<vector<int> > sPreviousObjectIDs_Gamma2; //1st dimension is cut-index
  • For each new event, clear the lists of previous object-ids
    • Put/call code like this at the beginning of your Process() function
//sNumCuts is the number of cuts you're placing in your TSelector
if((RunNumber != sPreviousRunNumber) || (EventNumber != sPreviousEventNumber))
{
  sPreviousObjectIDs_Beam.clear();
  sPreviousObjectIDs_PiMinus1.clear();
  sPreviousObjectIDs_PiMinus2.clear();
  sPreviousObjectIDs_PiPlus1.clear();
  sPreviousObjectIDs_PiPlus2.clear();
  sPreviousObjectIDs_Gamma1.clear();
  sPreviousObjectIDs_Gamma2.clear();
  sPreviousObjectIDs_Beam.resize(sNumCuts + 1);
  sPreviousObjectIDs_PiMinus1.resize(sNumCuts + 1);
  sPreviousObjectIDs_PiMinus2.resize(sNumCuts + 1);
  sPreviousObjectIDs_PiPlus1.resize(sNumCuts + 1);
  sPreviousObjectIDs_PiPlus2.resize(sNumCuts + 1);
  sPreviousObjectIDs_Gamma1.resize(sNumCuts + 1);
  sPreviousObjectIDs_Gamma2.resize(sNumCuts + 1);
}
  • For each cut you place when processing your event with your TSelector:
    • If it passes the cut: Increase the cut index
    • If it fails the cut: End the evaluation of the entry, and save the object-ids so that future entries can check against them.
size_t locCutIndex = 0;
if(Gamma1__Energy_BCAL < 0.5)
  return End_Entry(locCutIndex);
++locCutIndex;
 
Bool_t MySelector::End_Entry(size_t locCutIndex)
{
  sPreviousRunNumber = RunNumber;
  sPreviousEventNumber = EventNumber;
  for(size_t loc_i = 0; loc_i <= locCutIndex; ++loc_i)
  {
    sPreviousObjectIDs_Beam[loc_i].push_back(Beam__ObjectID);
    sPreviousObjectIDs_PiMinus1[loc_i].push_back(PiMinus1__ObjectID);
    sPreviousObjectIDs_PiMinus2[loc_i].push_back(PiMinus2__ObjectID);
    sPreviousObjectIDs_PiPlus1[loc_i].push_back(PiPlus1__ObjectID);
    sPreviousObjectIDs_PiPlus2[loc_i].push_back(PiPlus2__ObjectID);
    sPreviousObjectIDs_Gamma1[loc_i].push_back(Gamma1__ObjectID);
    sPreviousObjectIDs_Gamma2[loc_i].push_back(Gamma2__ObjectID);
  }
  return kTRUE;
}
  • Finally, when filling your histogram, first loop over the object-ids of any previous particle combinations to make sure that the result won't be duplicated:
bool locDuplicateFlag = false;
for(size_t loc_i = 0; loc_i < sPreviousObjectIDs_Gamma1[locCutIndex].size(); ++loc_i)
{
  if((sPreviousObjectIDs_Gamma1[locCutIndex][loc_i] == Gamma1__ObjectID) && (sPreviousObjectIDs_Gamma2[locCutIndex][loc_i] == Gamma2__ObjectID))
  {
    locDuplicateFlag = true;
    break;
  }
  if((sPreviousObjectIDs_Gamma1[locCutIndex][loc_i] == Gamma2__ObjectID) && (sPreviousObjectIDs_Gamma2[locCutIndex][loc_i] == Gamma1__ObjectID))
  {
    locDuplicateFlag = true; //g1 & g2 are switched but the invariant mass will be identical
    break;
  }
}
if(!locDuplicateFlag)
  MyPiZeroMassHist->Fill(locPiZeroMass);

Converting for AmpTools

  • To convert the TTree for use as input to AmpTools, use the sim-recon/src/programs/Utilities/tree_to_amptools/ program. Run it with no arguments to get usage instructions.

Retrieved from "https://halldweb1.jlab.org/wiki/index.php?title=Analysis_TTreeFormat&oldid=69948"