Algorithm Development Pt. 2

Analyzing tracks

Let's look at some tracks. In our execute() function, we'll use the selectedTracks() function (a member function belonging to our mother class xTRT::Algorithm). We're using it to retrieve, you guessed it, some selected tracks. What this function does is create a new container of xAOD::TrackParticle objects, where all of the tracks in the container passed the Tracks. selection that you define in your config file. You can go to the Config section for more detailed info, but here I'll just mention some track specific options. The default config we're using have options of the form: (the numbers might not be the same, but you can change them!)

### Cuts for the selectedTracks() container
Tracks.p: 5              ## momentum requirement in GeV
Tracks.pT: 5             ## transverse momentum requirement in GeV
Tracks.eta: 2.0          ## absolute value eta max
Tracks.nSi: 6            ## number of silicon hits required (pixel + SCT)
Tracks.nPix: 1           ## number of pixel hits required
Tracks.nTRT: 12          ## total number of TRT hits required
Tracks.nTRTprec: 1       ## total number of TRT precision hits required

Alright, so how would we get tracks that pass this selection? In just four lines, here's grabbing the selected tracks, looping over them, and grabbing the track occupancy (don't add this code yet, we'll add it in a sec.):

EL::StatusCode AnaProjectLooperAlg::execute() {
  // ...
  // ...

  auto tracks = selectedTracks();
  for ( const xAOD::TrackParticle* track : *tracks ) {
    float trkOccupancy = xTRT::trackOccupancy(track);
  }

  // ...
  // ...
}

trackOccupancy(...) is a helper function that is part of the TRTFramework library. You can find technical details in the doxygen API documentation here. You can find a list of helper functions that take either xAOD::TrackParticle, xAOD::Electron, or xAOD::Muon objects specifically here.

Let's make a tree that stores a few track properties and fills for every track as well as also some histograms. Go back to the header and define a few new private member variables:

private:
// ...

  TTree* m_trkTree;     //!
  float m_trkOccupancy; //!
  float m_trkpt;        //!
  int   m_trkNTRThits;  //!

And in our source:

EL::StatusCode AnaProjectLooperAlg::histInitialize() {
  // ..
  // ..
  // ..

  // new histograms
  create(TH1F("h_trkOccupancy","h_trkOccupancy",100,0,1));
  create(TH1F("h_trkpt",       "h_trkpt",       50,0,100));
  create(TH1F("h_trkNTRThits", "h_trkNTRThits", 50,0,50));

  // new tree output
  SETUP_OUTPUT_TREE(m_trkTree,"track_tree");
  m_trkTree->Branch("trkOccupancy",&m_trkOccupancy);
  m_trkTree->Branch("trkpt",       &m_trkpt);
  m_trkTree->Branch("trkNTRThits", &m_trkNTRThits);
  // let's also add the m_eventWeight variable we've already made to the event tree output
  m_trkTree->Branch("eventWeight",&m_eventWeight);

  return EL::StatusCode::SUCCESS;  
}

EL::StatusCode AnaProjectLooperAlg::execute() {
  // ...
  // ...
  // ...

  auto tracks = selectedTracks();
  // if you want all raw tracks from the sample, you can do:
  // auto tracks = trackContainer();
  for ( const xAOD::TrackParticle* track : *tracks ) {
    m_trkOccupancy = xTRT::trackOccupancy(track);
    m_trkpt = track->pt()*toGeV; // convert MeV to GeV
    m_trkNTRThits = xTRT::nTRT(track); // another helper function from mother
    m_trkTree->Fill(); // fill the tree
    // fill the histograms
    grab<TH1F>("h_trkOccupancy")->Fill(m_trkOccupancy,m_eventWeight);
    grab<TH1F>("h_trkpt")->Fill(m_trkpt,m_eventWeight);
    grab<TH1F>("h_trkNTRThits")->Fill(m_trkNTRThits,m_eventWeight);
  } // end loop over tracks

  // ..
  // ..

  return EL::StatusCode::SUCCESS;
}

Now go recompile, run the code, and find your new histograms and new tree + branches!

Analyzing hits

To analyze using hits, we need to access them. While looping over a track, we can grab the track measurement container like so: 1. create some pointers for the measurement on surface (msos) and the drift circle (driftCircle) 2. Use a predefined auxiliary accessor (xTRT::Acc::msosLink, which is an SG::AuxElement::ConstAccessor for the link to the measurement on surface) to check if the container of measurements is there, then loop over them. To read more about aux data go to that section. 3. for each track measurement, set the msos and driftCircle variables after making sure that their links are valid (also make sure we only have TRT hits).

The code block below shows the process.

At that point, you have access to the xTRT::MSOS (a typedef to xAOD::MeasurementOnSurface) object (we're calling our variable msos) and the xTRT::DriftCircle (a typedef to xAOD::TrackMeasurementValidation) object (we're calling our variable driftCircle).*

You can use the xTRT::getHitSummary(...) function to fill a struct with hit properties (instead of manually grabbing the aux data). See all the properties in the xTRT::HitSummary struct here. To learn more details about aux data, go to this chapter.

// inside execute()

// ...
bool isMCevent = isMC();
// ...

for ( const xAOD::TrackParticle* track : *tracks ) {
  // inside the loop over tracks

  // step 1 in above explanation
  const xTRT::MSOS* msos = nullptr;
  const xTRT::DriftCircle* driftCircle = nullptr;

  // step 2 in above explanatation
  if ( xTRT::Acc::msosLink.isAvailable(*track) ) {
    // loop described in step 3 above
    for ( const auto& trackMeasurementLink : xTRT::Acc::msosLink(*track) ) {
        // step 3 continued, make sure the link to the track measurement is valid
        if ( trackMeasurementLink.isValid() ) {
        // get the track measurement by dereferencing the link
        msos = *trackMeasurementLink;
        if ( msos->detType() != 3 ) continue; // makes sure we get TRT hits only.
        // make sure the drift circle (the track measurement validation) link is valid
        if ( !(msos->trackMeasurementValidationLink().isValid()) ) continue;
        // get the drift circle by dereferencing the link
        driftCircle = *(msos->trackMeasurementValidationLink());

        // using the getHitSummary function:
        auto  hitSummary = xTRT::getHitSummary(track,msos,driftCircle,isMCevent);
        bool  hit_is_HT  = hitSummary.HTMB;
        float hit_tot    = hitSummary.ToT;
        /// ...

      }
    }
  }
  else {
    ANA_MSG_WARNING("msosLink was not available!");
  }

* We define these typedef's because it's more intuitive to use MSOS and DriftCircle. The actual classes are used throughout the inner detector software, that's why they have a generic name.