Difference between revisions of "Analysis TTreeFormat"

From GlueXWiki
Jump to: navigation, search
(Accessing TTree Data)
(Accessing TTree Data)
Line 276: Line 276:
 
<br style="clear:both;"/>
 
<br style="clear:both;"/>
 
* <span style="color:#0000FF">TSelector</span> (histogram b1pi mass distributions):
 
* <span style="color:#0000FF">TSelector</span> (histogram b1pi mass distributions):
** [http://root.cern.ch/drupal/content/developing-tselector Link: TSelector Documentation]
+
** [http://root.cern.ch/drupal/content/proof Documentation Link: PROOF]
 +
** [http://root.cern.ch/drupal/content/proof-multicore-desktop-laptop-proof-lite Documentation Link: PROOF-Lite]
 +
** [http://root.cern.ch/drupal/content/developing-tselector Documentation Link: TSelector]
 +
** [http://root.cern.ch/drupal/content/progress-dialog Misc. Documentation Link: Disable Graphics]
 +
** [http://root.cern.ch/drupal/content/processing-proof Documentation Link: Full TSelector Example (with PROOF-Lite)]
 +
** [http://root.cern.ch/drupal/content/basic-processing Documentation Link: Process Examples]  
 +
 
 
<syntaxhighlight>
 
<syntaxhighlight>
 
GetEntry(entry);
 
GetEntry(entry);

Revision as of 22:51, 11 September 2013

Why a standard TTree format?

  • For any reaction can streamline (and provide a best-practices implementation/examples of): analysis utilities, BDT setup/input, amplitude analysis setup/input(?)
  • Makes it easy for users to keep everything organized, especially handling of the combinatoric background.
  • Format is designed to be one-size-fits-all, but is still extensible (customizable).

TTree Format - Overview

ROOT Version

  • ROOT >= v5.32 will be required for building DANA
    • Much easier to code TTree creation: TClonesArray::ConstructedAt() can be used instead of calls to placement-new
    • Just about any version of ROOT is fine for reading the data

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 particle combination.
  • One TBranch per variable. (e.g. "RunNumber", "PiMinus__ObjectID", etc.)
  • Thrown particles and unused particle hypotheses will be stored in arrays/TClonesArray's: one array index per particle.
    • e.g.: "Unused__ChiSq_Tracking" (Double_t[]), "Thrown__P4_Thrown" (TClonesArray(TLorentzVector))
    • 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

  • 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 - Detail

Branch Name Prefixes

Example Reaction (b1pi):

  • γ p → X(2000), (p)
    • X(2000) → b1(1235)+, π-
      • b1(1235)+ → ω, π+
        • ω → π+, π-, π0
          • π0 → γ γ

Branch Names:

  • Beam: "BeamPhoton"
  • Detected Reaction Particles: "PiMinus1," "PiPlus1," "PiPlus2," "PiMinus2," "Photon1," "Photon2"
  • Decaying Reaction Particles: "DecayingPi0"
    • No branches are created for resonances.
  • Missing: "MissingProton"
  • Thrown Beam: "ThrownBeam"
  • Thrown Products: "Thrown"
    • Array entries do NOT correspond to particle combos: just a different particle in each array index
    • Note: particles with PID = Unknown, particles decaying from final-state particles, and orphan particles are not saved.
  • Unused Particle Hypothesis Name: "Unused"
    • Array entries do NOT correspond to particle combos: just a different particle hypothesis in each array index

Non-Particle Data

// EVENT DATA
"RunNumber": UInt_t
"EventNumber": UInt_t
 
// # PARTICLES //these are the array sizes for many of the other branches
"NumThrown": UInt_t
"NumUnused": UInt_t
 
// RF
"RFTime_Thrown": Double_t
"RFTime_Measured": Double_t
"RFTime_KinFit": Double_t //only if kinematic fit performed
 
// KINEMATIC FIT
"ChiSq_KinFit": Double_t //only if kinematic fit performed
"NDF_KinFit": UInt_t //only if kinematic fit performed // = 0 if kinematic fit doesn't converge
 
//MONTE CARLO INFORMATION
"MCWeight": Double_t
"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()

Particle Data : Thrown Beam Particles

  • Full Branch Name Example: "ThrownBeam__PID"
//IDENTIFIER
"PID": Int_t //PDG ID value
 
//KINEMATICS: THROWN  //At the production vertex 
"X4_Thrown": TLorentzVector
"P4_Thrown": TLorentzVector

Particle Data : Thrown Products

  • Note: particles with PID (Particle_t) = Unknown, particles decaying from final-state particles, and orphan particles are not saved.
  • Full Branch Name Example: "Thrown__ParentID"
//IDENTIFIERS / MATCHING
"ParentID": 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
"MatchID": Int_t["NumThrown"] //the "ObjectID" of the unused and/or detected reaction particle it is matched with (-1 for no match)
 
//KINEMATICS: THROWN  //At the production vertex 
"X4_Thrown": TClonesArray(TLorentzVector["NumThrown"])
"P4_Thrown": TClonesArray(TLorentzVector["NumThrown"])

Particle Data : Beam Reaction Particles

  • Full Branch Name Example: "Beam__ObjectID"
//IDENTIFIER
"ObjectID": Int_t //each beam particle has its own #
 
//KINEMATICS: MEASURED //At the target center
"X4_Measured": TLorentzVector
"P4_Measured": TLorentzVector
 
//KINEMATICS: KINFIT //At the interaction vertex //only present if kinfit performed
"X4_KinFit": TLorentzVector
"P4_KinFit": TLorentzVector

Particle Data : Detected Reaction Particles

  • Full Branch Name Example: "PiMinus1__ObjectID"
//IDENTIFIER / MATCHING
"ObjectID": Int_t //each physical particle has its own # (to keep track of different pid hypotheses for the same particle)
"MatchID": Int_t //the array index of the thrown particle it is matched with (-1 for no match)
 
//KINEMATICS: MEASURED //At the production vertex
"X4_Measured": TLorentzVector //t is the measured value in TOF/BCAL/FCAL projected back to Position_Measured
"P4_Measured": TLorentzVector
 
//KINEMATICS: KINFIT //At the production vertex //only present if kinfit performed
"X4_KinFit": TLorentzVector
"P4_KinFit": TLorentzVector
 
// PID QUALITY:
"NDF_Tracking": UInt_t //for charged only
"ChiSq_Tracking": Double_t //for charged only
"NDF_Timing": UInt_t
"ChiSq_Timing_Measured": Double_t //using measured data
"ChiSq_Timing_KinFit": Double_t //using kinematic fit data //only present if kinfit performed
"NDF_DCdEdx": UInt_t //for charged only
"ChiSq_DCdEdx": Double_t //for charged only
 
// DEPOSITED ENERGY:
"dEdx_CDC": Double_t //for charged only
"dEdx_FDC": Double_t //for charged only
"dEdx_TOF": Double_t //for charged only
"dEdx_ST": Double_t //for charged only
"Energy_BCAL": Double_t
"Energy_FCAL": Double_t

Particle Data : Decaying Reaction Particles

  • Excludes resonances
  • Full Branch Name Example: "DecayingPi0__P4"
//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": TLorentzVector //only present if kinfit performed

Particle Data : Missing Reaction Particles

  • Full Branch Name Example: "MissingProton__P4"
//KINEMATICS: //At the decay vertex 
"X4": TLorentzVector //kinematic fit result if kinfit performed, else reconstructed from detected particles
"P4": TLorentzVector //only present if kinfit performed

Particle Data : Unused Hypotheses

  • Full Branch Name Example: "Unused__ObjectID"
//IDENTIFIERS / MATCHING
"ObjectID": Int_t["NumUnused"] //each physical particle has its own # (to keep track of different pid hypotheses for the same particle)
"PID": Int_t["NumUnused"] //PDG ID value
"MatchID": Int_t //the array index of the thrown particle it is matched with (-1 for no match)
 
//KINEMATICS: MEASURED  //At the production vertex 
"P4_Measured": TClonesArray(TLorentzVector["NumUnused"])
"X4_Measured": TClonesArray(TLorentzVector["NumUnused"]) //t is the measured value in TOF/BCAL/FCAL projected back to Position_Measured
 
// PID QUALITY:
"NDF_Tracking": UInt_t["NumUnused"] //for charged only
"ChiSq_Tracking": Double_t["NumUnused"] //for charged only
"NDF_Timing": UInt_t["NumUnused"]
"ChiSq_Timing": Double_t["NumUnused"]
"NDF_DCdEdx": UInt_t["NumUnused"] //for charged only
"ChiSq_DCdEdx": Double_t["NumUnused"] //for charged only
 
// DEPOSITED ENERGY:
"dEdx_CDC": Double_t["NumUnused"] //for charged only
"dEdx_FDC": Double_t["NumUnused"] //for charged only
"dEdx_TOF": Double_t["NumUnused"] //for charged only
"dEdx_ST": Double_t["NumUnused"] //for charged only
"Energy_BCAL": Double_t["NumUnused"]
"Energy_FCAL": Double_t["NumUnused"]

Event-Independent DReaction Information

  • 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.
#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

  • This is basically something that you just have to do manually.
    • Note: this actually isn't thread safe: another thread can modify the value to be written to the tree in between the last unlock below and the tree fill call...
//Create the trees
vector<const DEventWriterROOT*> locEventWriterROOTVector;
locEventLoop->Get(locEventWriterROOTVector); //creates the TTrees for all DReactions upon first call
 
//Create the branch
Double_t* locMyVariable = new Double_t;
string locReactionName = locReaction->Get_ReactionName();
string locTreeName = locReactionName + string("_Tree");
japp->RootWriteLock(); //always acquire a lock before accessing the global ROOT file
{
  //get the tree
  gDirectory->cd("/");
  gDirectory->cd(locReactionName.c_str());
  TTree* locTree = (TTree*)gDirectory->Get(locTreeName.c_str());
  locTree->Branch("MyBranchName", locMyVariable, "D");
}
japp->RootUnLock();
 
//Set the data in the branch and fill
japp->RootWriteLock(); //lock: any thread has access to the tbranch
{
  *locMyVariable = 3.0; //or whatever you want to save
}
japp->RootUnLock();
locEventWriterROOTVector[0]->Fill_Trees(locEventLoop);

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
//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
//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.