---+ %NOTOC% xAOD Analysis in !EventLoop

Let's setup the ATLAS software environment:

<verbatim style="background: #e0ebf6;">

setupATLAS

</verbatim>

This alias is already defined for you when you log in to lxplus. After you type it you will see a list of commands you can type to setup various ATLAS computing tools.

<blockquote>

If you are not working on lxplus: you will need to define these variables: <br/>

<verbatim>

export ATLAS_LOCAL_ROOT_BASE=/cvmfs/atlas.cern.ch/repo/ATLASLocalRootBase

alias setupATLAS='source ${ATLAS_LOCAL_ROOT_BASE}/user/atlasLocalSetup.sh'

setupATLAS

</verbatim>

</blockquote>

We will use the ATLAS Analysis Release (new for Run 2) to setup our environment and other RootCore packages with correct package tags (such as those recommended by the CP groups). You can use the Analysis Release if you have access to cvmfs (lxplus, local clusters, machines, and the Grid).

<blockquote>

It is also possible to checkout all packages in a given Analysis Release and manually compile them on your working machine. This alternative setup is described below. But for this tutorial we recommend you just use the release from cvmfs on lxplus.

</blockquote>

* =finalize()=

%SYNTAX{ syntax="cpp" }%

float muPtCone20 = 0.; // your variable that will be filled after calling the isolation function

(*muon_itr)->isolation(muPtCone20, xAOD::Iso::ptcone20); // second arg is an enum defined in xAODPrimitives/IsolationType.h

%ENDSYNTAX%

Alternatively you can access that same variable by simply doing:

%SYNTAX{ syntax="cpp" }%

(*muon_itr)->auxdata< float >("ptcone20");

%ENDSYNTAX%

For the muons you can find the complete list of accessors in the xAOD Muon class (version 1)

---+ 5. Creating and running our steering macro

To actually run this EventLoop algorithm we need some steering code. This can be a root macro in either C++ or python or some compiled C++ code. For this tutorial we will create a C++ macro.

We will use another ASG tool called SampleHandler which is a nice tool that allows for easy sample management. In this example we will create and configure a SampleHandler object. We will specify the path to the main directory, under which there could be several subdirectories (typically representing datasets) and within those the individual input files. Here we will tell SampleHandler we are only interested in one input xAOD file (specified by the exact name, but wildcards are accepted to use several specific inputs). More information and options for using SampleHandler to 'find' your data is found on the dedicated SampleHandler wiki.

First let's create a Run directory for it. From your main working directory execute:

<verbatim style="background: #e0ebf6;">

mkdir Run

cd Run

</verbatim>

And in that directory we will create a new file called ATestRun.cxx. Fill this new file, ATestRun.cxx, with the following:

%SYNTAX{ syntax="cpp" }%

void ATestRun (const std::string& submitDir)

{

//===========================================

// FOR ROOT6 WE DO NOT PUT THIS LINE

// (ROOT6 uses Cling instead of CINT)

// Load the libraries for all packages

// gROOT->Macro("$ROOTCOREDIR/scripts/load_packages.C");

// Instead on command line do:

// > root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

// The above works for ROOT6 and ROOT5

//==========================================

// Set up the job for xAOD access:

xAOD::Init().ignore();

// create a new sample handler to describe the data files we use

SH::SampleHandler sh;

// scan for datasets in the given directory

// this works if you are on lxplus, otherwise you'd want to copy over files

// to your local machine and use a local path. if you do so, make sure

// that you copy all subdirectories and point this to the directory

// containing all the files, not the subdirectories.

// use SampleHandler to scan all of the subdirectories of a directory for particular MC single file:

job.algsAdd (alg);

// make the driver we want to use:

// this one works by running the algorithm directly:

EL::DirectDriver driver;

// we can use other drivers to run things on the Grid, with PROOF, etc.

// process the job using the driver

driver.submit (job, submitDir);

}

%ENDSYNTAX%

Read over the comments carefully to understand what is happening. Notice that we will only run over the first 500 events (for testing purposes). Obviously if you were doing a real analysis you would want to remove that statement to run over all events in a sample.

Ok, now the big moment has come. Within your Run directory execute your ATestRun.cxx macro with root:

<verbatim style="background: #e0ebf6;">

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

You can quit ROOT by typing ".q" at the ROOT prompt.

Note that submitDir is the directory where the output of your job is stored. If you want to run again, you either have to remove that directory or pass a different name into ATestRun.cxx.

<blockquote>

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

</verbatim>

---+++ The Good Runs List

The Good Runs List (GRL) is an xml file that selects good luminosity blocks from within the data runs (spanning 1-2 minutes of data-taking). Luminosity blocks which are not listed in this xml file are considered bad, and should not be used in your analysis.

Now let's move to our source code MyAnalysis/Root/MyxAODAnalysis.cxx, and initialize the tool for usage in our initialize() method:

%SYNTAX{ syntax="cpp" }%

// GRL

m_grl = new GoodRunsListSelectionTool("GoodRunsListSelectionTool");

}

} // end if not MC

%ENDSYNTAX%

And last but not least, remember to clean up the memory in finalize():

%SYNTAX{ syntax="cpp" }%

// GRL

if (m_grl) {

delete m_grl;

m_grl = 0;

}

%ENDSYNTAX%

Finally let's compile! Note we have to run rc find_packages again, as we updated the package dependencies in our Makefile.

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

</verbatim>

Choices of GRLs are described on this page: Good Run Lists for Analysis

To add: Determine the integrated luminosity of our data sample <br/>

---+++ Removing additional detector imperfections

The GRL helps remove detector problems that affect entire luminosity blocks (1-2 minutes of data-taking). Note that even one luminosity block can be thousands of events. To avoid throwing away perfectly fine events there are some additional event-level detector flags, to help single events that are plagued with some detector problem. For the 2015 data-taking period there are three flags to remove problematic events:

* due to the liquid argon system

* due to the tile calorimeter system

* due to the SCT inner detector system

* due to incomplete events (event information missing after TTC restarts)

Here we will show you how to remove these problematic events. Here we assume the EventInfo EDM classes have already been defined and used from instructions at the top of this page.

We can count the number of 'clean' events, by defining a member variable in our header file MyAnalysis/MyAnalysis/MyxAODAnalysis.h:

%SYNTAX{ syntax="cpp" }%

int m_numCleanEvents; //!

%ENDSYNTAX%

Now in our source code let's initialize this event count to zero, so in the initialize() method in MyAnalysis/Root/MyxAODAnalysis.cxx simply add this line:

%SYNTAX{ syntax="cpp" }%

m_numCleanEvents = 0;

%ENDSYNTAX%

Now for the cleaning! This is something done event-by-event (think execute() method!), and only for data:

%SYNTAX{ syntax="cpp" }%

//------------------------------------------------------------

// Apply event cleaning to remove events due to

// problematic regions of the detector

// or incomplete events.

// Apply to data.

//------------------------------------------------------------

// reject event if:

if(isMC){

if( (eventInfo->errorState(xAOD::EventInfo::LAr)==xAOD::EventInfo::Error ) || (eventInfo->errorState(xAOD::EventInfo::Tile)==xAOD::EventInfo::Error ) || (eventInfo->errorState(xAOD::EventInfo::SCT)==xAOD::EventInfo::Error ) || (eventInfo->isEventFlagBitSet(xAOD::EventInfo::Core, 18) ) )

{

return EL::StatusCode::SUCCESS; // go to the next event

} // end if event flags check

} // end if the event is data

m_numCleanEvents++;

%ENDSYNTAX%

You can print out the final number of clean events from the finalize() method:

%SYNTAX{ syntax="cpp" }%

Info("finalize()", "Number of clean events = %i", m_numCleanEvents);

%ENDSYNTAX%

This will do nothing for MC files, but could actually remove some events from the 2015 data samples (you could try using the data xAOD file as input in your steering macro).

You can compile your package as before (we do not need to run rc find_packages again as we did not update the package dependencies), and you can run your job in ROOT from the Run/ directory as before. From your working directory:

<verbatim style="background: #e0ebf6;">

rc compile

cd Run/

rm -rf submitDir/

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

---++ Jets

Let's create a basic loop over a jet container, for the AntiKt4EMTopoJets jet collection:

First, let's add the jet xAOD EDM package to MyAnalysis/cmt/Makefile.RootCore:

%SYNTAX{ syntax="cpp" }%

PACKAGE_DEP = EventLoop [...] xAODJet

%ENDSYNTAX%

where [...] is any other package dependencies you may already have included.

Now in our source code MyAnalysis/Root/MyxAODAnalysis.cxx, add the jet xAOD EDM container class to include statement near the top:

%SYNTAX{ syntax="cpp" }%

#include "xAODJet/JetContainer.h"

%ENDSYNTAX%

Still in the source code, in execute() (run for every event) retrieve the jet container and loop over all jets for this event:

%SYNTAX{ syntax="cpp" }%

// get jet container of interest

%ENDSYNTAX%

If you are not sure of the container type and/or container name see above.

You can do the usual RootCore thing to update/find the (new) package dependencies and compile:

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

</verbatim>

---+++ Jet cleaning tool

Let's add the jet cleaning tool to our code. This jet selector tool applies the jet cleaning described on this page: How to clean jets 2015. It is a cleaning that you could apply to jets in both data and MC (see the dedicated wiki page for the recommended use).

Below describe the steps necessary to including this jet cleaning in our analysis code:

First, add the package dependency to MyAnalysis/cmt/Makefile.RootCore:

%SYNTAX{ syntax="cpp" }%

PACKAGE_DEP = EventLoop [...] JetSelectorTools

%ENDSYNTAX%

where [...] is any other package dependencies you may already have included.

Second, add the class dependency to our header file, at the top near where all the other "includes" are, we will add:

%SYNTAX{ syntax="cpp" }%

#include "JetSelectorTools/JetCleaningTool.h"

%ENDSYNTAX%

Still in the header file, let's create an instance of the tool. In the class definition (so somewhere within class MyxAODAnalysis : public EL::Algorithm{...}) :

%SYNTAX{ syntax="cpp" }%

JetCleaningTool *m_jetCleaning; //!

%ENDSYNTAX%

Now moving to our source code,=MyAnalysis/Root/MyxAODAnalysis.cxx= ,we need to initialize this tool, in the initialize() method add these lines:

%SYNTAX{ syntax="cpp" }%

// initialize and configure the jet cleaning tool

m_jetCleaning = new JetCleaningTool("JetCleaning");

m_jetCleaning->msg().setLevel( MSG::DEBUG );

ANA_CHECK(m_jetCleaning->setProperty( "CutLevel", "LooseBad"));

ANA_CHECK(m_jetCleaning->setProperty("DoUgly", false));

ANA_CHECK(m_jetCleaning->initialize());

%ENDSYNTAX%

Then, in finalize(), don't forget to delete the tool from memory.

%SYNTAX{ syntax="cpp" }%

if( m_jetCleaning ) {

delete m_jetCleaning;

m_jetCleaning = 0;

}

%ENDSYNTAX%

In the above we have set a DEBUG level of text information to be printed out, and we have set the cut level to "MediumBad", there are several other options; the choice however is usually analysis dependent.

%ENDCODE%

After the jet loop you could print out the number of good jets.

<blockquote>

Apart from the setup and configuration, here is a quick comparison of how this tool was called for jets using D3PDs and how it is now called using xAODs:

Ta-da! Isn't that just SO much nicer :-)

</blockquote>

---+++ Jet energy resolution

This tool takes jets as inputs and smears the MC such that the energy resolution matches that of data (note in 2015 the data and MC agreed pretty well, so this smearing was not necessary). The tool also provides the uncertainty associated with the energy resolution. More official information about this tool (including configuration options) can be found on the following jet pages:

<!--Jet energy resolution provider 2012-->

* JetEnergyResolutionXAODTools

* JetResolution2015Prerecom

Here we will show you how to implement the tool. Like all tools, there are several steps:

double uncert = m_JERTool->getUncertainty(*jet_itr); // you can provide a second argument specify which nuisance parameter, the default is all

Info("execute()", "jet mcRes = %f , uncert = %f", mcRes, uncert);

} // end if MC

...

// in finalize, delete tool:

if(m_JERTool){

delete m_JERTool;

m_JERTool = 0;

}

%ENDSYNTAX%

We've done a lot of stuff - let's recompile the package and test that we didn't break anything:

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

cd Run/

rm -rf submitDir/

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

---++ Muons

In this introductory part we will show you how to make a basic loop over muon objects. This is very much the same thing we did for jets already.

Add the muon xAOD EDM package to MyAnalysis/cmt/Makefile.RootCore:

%SYNTAX{ syntax="cpp" }%

PACKAGE_DEP = EventLoop [...] xAODMuon

%ENDSYNTAX%

where [...] is any other package dependencies you may already have included.

Now in our source code MyAnalysis/Root/MyxAODAnalysis.cxx, add the muon xAOD EDM container class to include statement near the top:

%SYNTAX{ syntax="cpp" }%

#include "xAODMuon/MuonContainer.h"

%ENDSYNTAX%

Still in the source code, in execute() (run for every event) retrieve the muon container and loop over all muons for this event:

%SYNTAX{ syntax="cpp" }%

//------------

// MUONS

//------------

// get muon container of interest

%ENDSYNTAX%

If you are not sure of the container type and/or container name see above.

You can do the usual RootCore thing to update/find the (new) package dependencies and compile:

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

</verbatim>

More information about using muons for DC14 can be found here: <br/>

MCP: Analysis Guidelines for DC14

---+++ Muon calibration and smearing tool, and associated systematics

Now let's use the MuonCalibrationAndSmearingTool found in the MuonMomentumCorrections package which is used to correct the MC to look like the data.

As with all tools we have to do the basics:

xAOD::MuonContainer::iterator muonSC_itr = (muons_shallowCopy.first)->begin();

xAOD::MuonContainer::iterator muonSC_end = (muons_shallowCopy.first)->end();

for( ; muonSC_itr = muonSC_end; ++muonSC_itr ) {

if(m_muonCalibrationAndSmearingTool->applyCorrection(**muonSC_itr) = CP::CorrectionCode::Error){ // apply correction and check return code

// Can have CorrectionCode values of Ok, OutOfValidityRange, or Error. Here only checking for Error.

// If OutOfValidityRange is returned no modification is made and the original muon values are taken.

Error("execute()", "MuonCalibrationAndSmearingTool returns Error CorrectionCode");

}

Info("execute()", " corrected muon pt = %.2f GeV", ((*muonSC_itr)->pt() * 0.001));

} // end for loop over shallow copied muons

delete muons_shallowCopy.first;

delete muons_shallowCopy.second;

%ENDSYNTAX%

%SYNTAX{ syntax="cpp" }%

// header for systematics:

#include "PATInterfaces/SystematicRegistry.h"

...

// list of systematics

std::vector<CP::SystematicSet> m_sysList; //!

%ENDSYNTAX%

// apply recommended systematic for muonCalibrationAndSmearingTool

if( m_muonCalibrationAndSmearingTool->applySystematicVariation( *sysListItr ) = CP::SystematicCode::Ok ) {

Error("execute()", "Cannot configure muon calibration tool for systematic" );

continue; // go to next systematic

} // end check that systematic applied ok

// create a shallow copy of the muons container

xAOD::MuonContainer::iterator muonSC_itr = (muons_shallowCopy.first)->begin();

xAOD::MuonContainer::iterator muonSC_end = (muons_shallowCopy.first)->end();

for( ; muonSC_itr = muonSC_end; ++muonSC_itr ) {

m_muonCalibrationAndSmearingTool->applyCorrection(*muonSC_itr);

Info("execute()", "corrected muon pt = %.2f GeV", ((*muonSC_itr)->pt() * 0.001));

} // end for loop over shallow copied muons

delete muons_shallowCopy.first;

delete muons_shallowCopy.second;

} // end for loop over systematics

%ENDSYNTAX%

Now you can recompile to see this new tool in action:

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

cd Run/

rm -rf submitDir/

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

---++ Taus

More information about looping over taus and accessing tau xAOD quantities can be found on this page: <br/>

Tau xAOD Instructions

Information about tau-related tools for:

TauSelectionTool: generic tool to apply a set of requirements on tau candidates

* TauSmearingTool: currently supports tau energy corrections

* TauEfficiencyCorrectionsTool: currently provides jet identification scale factors and the associated uncertainties

can be read in this very nice documentation included in the package: <br/>

https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/TauID/TauAnalysisTools/tags/TauAnalysisTools-00-00-04/README.rst <br/>

(Note I have linked to TauAnalysisTools-00-00-04 as this is the tagged version in AnalysisBase 2.0.6., but this can/will change with release)

---++ MC Truth

Lots of goodies about how to extract MC truth information in this talk: Monte Carlo truth in the xAOD

---++ Trigger

Here we will show you how to do some simple things with the trigger (namely if the trigger passed/failed and the final trigger prescale factor). First let me mention that the trigger is a complicated beast... It relies on meta-data, which is "data about the data", like which trigger menu was used which contains info about the thresholds, chains, and prescales. And it relies on event-by-event level information, like which triggers passed/failed and the associated trigger objects that fired these triggers (like the object that was identified as an electron at the trigger level that fired a particular electron trigger). Because of this level of complication you can't just browse the trigger decision in a TBrowser in ROOT.

The primary tool you must use to interact with the trigger is called the TrigDecisionTool. Note that if you want to use the TrigDecisionTool in EventLoop (which is what we are doing here in this tutorial), you must also make use of another tool: the TrigConf::xAODConfigTool. Both of which must be declared on the "heap" (meaning a member variable declared in your header file). But we'll lead you through all that below!

Let's do the basics :

%ENDSYNTAX%

Now inside the source file MyAnalysis/Root/MyxAODAnalysis.cxx, we need to book the histogram and add it to the output. That happens in the histInitialize function:

%SYNTAX{ syntax="cpp" }%

EL::StatusCode MyxAODAnalysis :: histInitialize ()

{

// Here you do everything that needs to be done at the very

// beginning on each worker node, e.g. create histograms and output

// trees. This method gets called before any input files are

// connected.

h_jetPt = new TH1F("h_jetPt", "h_jetPt", 100, 0, 500); // jet pt [GeV]

wk()->addOutput (h_jetPt);

return EL::StatusCode::SUCCESS;

}

%ENDSYNTAX%

This method is called before processing any events. Note that the wk()->addOutput call is a mechanism EventLoop uses for delivering the results of an algorithm to the outside world. When running in PROOF, ROOT will merge all of the objects in this list. In principle you can place any object that inherits from TObject in there, but if you put anything other than histograms there you need to take special precautions for the merging step.

Now that we have the histogram, we also need to fill it. For that we can use the execute() method (called for every event). You probably have a loop over

the AntiKt4EMTopoJets collection (as described above), in that jet loop let's fill this histogram, so the loop looks something like:

%SYNTAX{ syntax="cpp" }%

// loop over the jets in the container

xAOD::JetContainer::const_iterator jet_itr = jets->begin();

xAOD::JetContainer::const_iterator jet_end = jets->end();

for( ; jet_itr = jet_end; ++jet_itr ) {

Info("execute()", " jet pt = %.2f GeV", ((*jet_itr)->pt() * 0.001)); // just to print out something

h_jetPt->Fill( ( (*jet_itr)->pt()) * 0.001); // GeV

} // end for loop over jets

%ENDSYNTAX%

Now we should be able to build our package with the full algorithm:

<verbatim style="background: #e0ebf6;">

rc compile

</verbatim>

And rerun from the Run/ directory:

<verbatim style="background: #e0ebf6;">

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

%ENDSYNTAX%

Notice outputName does not have a //! beside it, that's because we will defined it at run time (in the macro). And before closing this header file we need to include the TTree header file near the top of our code:

%SYNTAX{ syntax="cpp" }%

#include <TTree.h>

%ENDSYNTAX%

Now in our source code, MyAnalysis/Root/MyxAODAnalysis.cxx, let's do the actual implementation. First in the histInitialize() method, add these lines:

%SYNTAX{ syntax="cpp" }%

// get the output file, create a new TTree and connect it to that output

// define what braches will go in that tree

TFile *outputFile = wk()->getOutputFile (outputName);

tree = new TTree ("tree", "tree");

tree->SetDirectory (outputFile);

tree->Branch("EventNumber", &EventNumber);

%ENDSYNTAX%

Now in the execute() method let's actually fill the branch. You should do this somewhere after you have retrieved the EventInfo object:

%SYNTAX{ syntax="cpp" }%

// fill the branches of our trees

EventNumber = eventInfo->eventNumber();

%ENDSYNTAX%

And then somewhere near the end of the execute() method be sure to fill the tree:

%SYNTAX{ syntax="cpp" }%

tree->Fill();

%ENDSYNTAX%

Finally near the top of the source code add the include statement for the TFile object we are using:

%SYNTAX{ syntax="cpp" }%

#include <TFile.h>

%ENDSYNTAX%

To MyAnalysis/cmt/Makefile.RootCore add the dependency for EventLoopAlgs which where the NtupleSvc lives:

%SYNTAX{ syntax="cpp" }%

PACKAGE_DEP = EventLoop [...] EventLoopAlgs

%ENDSYNTAX%

After these modifications we should re-compile our package and make sure it works:

<verbatim style="background: #e0ebf6;">

rc find_packages

rc compile

</verbatim>

Moving now to our Run/ directory and our macro ATestRun.cxx, after you have defined the job but before you have added your algorithm to the job, add these lines which will create an output and an ntuple in that output stream:

%SYNTAX{ syntax="cpp" }%

// define an output and an ntuple associated to that output

EL::OutputStream output ("myOutput");

job.outputAdd (output);

EL::NTupleSvc *ntuple = new EL::NTupleSvc ("myOutput");

job.algsAdd (ntuple);

%ENDSYNTAX%

After this line job.algsAdd (alg); add the following line, which is basically letting your algorithm know the name of the output file:

%SYNTAX{ syntax="cpp" }%

alg->outputName = "myOutput"; // give the name of the output to our algorithm

%ENDSYNTAX%

(If you are using the compiled macro you will need to add two include statements for the NtupleSvc and Output ntuple, #include <EventLoopAlgs/NTupleSvc.h> and #include <EventLoop/OutputStream.h>.)

And that's it. You can rerun your macro again:

<verbatim style="background: #e0ebf6;">

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

Your new ntuple will appear in the directory submitDir/data-myOutput/.

---+ 9. Playing with xAOD containers

Here we will show you how to modify or create new xAOD containers/objects, and then how to write these to another xAOD in ROOT.

Note that if you do it this way you will not be able to read the output xAOD back into Athena. We will make use of some methods in the xAODRootAccess package to create our new xAOD.

%ENDSYNTAX%

---++ Copying full container(s) to a new xAOD

Here we will show you how to copy the contents of a full container, un-modified, for every event. We assume you have followed the instructions above to define a new output xAOD in EventLoop.

We will create this copy in the event loop, so in MyAnalysis/Root/MyxAODAnalysis.cxx in the execute() method add the following line to copy the full container for AntiKt4EMTopoJets:

%SYNTAX{ syntax="cpp" }%

// copy full container(s) to new xAOD

// without modifying the contents of it:

ANA_CHECK(event->copy("AntiKt4EMTopoJets"));

%ENDSYNTAX%

%ENDSYNTAX%

Compile like usual and test your code:

<verbatim style="background: #e0ebf6;">

rc compile

cd Run

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

(Remember that if you already have a submitDir directory the job will crash.)

If you have followed the instructions above you will find your output xAOD in submitDir/data-outputLabel/.

Note that you can only copy xAOD objects and/or containers. You can determine if the container is of xAOD type by running checkxAOD.py on the xAOD and seeing which objects are of type xAOD::x (where x is the container of interest).

---++ Creating modified container(s)

Sometimes you don't want to just blindly copy the entire contents of a container for every event into a new xAOD, sometimes you want to modify the (content of the) containers themselves. You could imagine you might want to do some of the following:

...

}

%ENDSYNTAX%

Also make sure you have updated your package dependencies in MyAnalysis/cmt/Makefile.RootCore to include:

%SYNTAX{ syntax="cpp" }%

PACKAGE_DEP = EventLoop [...] xAODJet xAODCore

%ENDSYNTAX%

where [...] represents any other package dependencies you may already have included.

Now to write this new container to a new xAOD output, we assume you have first setup the new xAOD output in EventLoop as we just described above:

Creating a new xAOD output in EventLoop

And then, still in execute() in the source code, after you have created and filled the new goodJet collection above, add these lines:

%SYNTAX{ syntax="cpp" }%

// Record the objects into the output xAOD:

ANA_CHECK(event->record( goodJets, "GoodJets" ));

ANA_CHECK(event->record( goodJetsAux, "GoodJetsAux." ));

%ENDSYNTAX%

Of course all of this needs to be put before the event->fill(); call somewhere. Also take note that you can call the interface container whatever you want (within reason, it should be an alphanumeric name), but the accompanying auxiliary store object must always have the same name, postfixed by "Aux.". So, for any "Bla" interface container you would record a "BlaAux." auxiliary container.

Notice that we create the output objects on the heap. (With new.) This mimics the behaviour of StoreGate. Just like StoreGate, TEvent also takes ownership of the objects recorded in it, so you must not delete the containers that you successfully recorded into TEvent.

You can compile and test your code:

<verbatim style="background: #e0ebf6;">

rc compile

cd Run

root -l '$ROOTCOREDIR/scripts/load_packages.C' 'ATestRun.cxx ("submitDir")'

</verbatim>

The new output xAOD in submitDir/data-outputLabel/ will now contain our new "GoodJets" container. You can check that it has fewer jets than the original AntiKt4EMTopoJets from which we deep copied.

---+++ Shallow copy

Another way of making copies of containers is called the shallow copy. If you want to apply some sort of modification/calibration to all objects in an input container, but you don't want to perform an object selection at the same time, then the best idea is to make a "shallow copy" of the input container. This is done with the help of xAOD::shallowCopyContainer. Note that you can not add/remove objects to/from a shallow copy container (but you can decorate it with new variables associated to each object - see decorations in the next section). However it's absolutely possible to make deep copies of some selected shallow copies later on in your code, and put these deep copies into a container of selected objects.

This creates a copy which only overrides default values with new ones that you set, while keeping all other unaffected quantities unchanged from the original. The shallowCopyContainer will return a pair of xAOD objects (one for the interface and one for the auxiliary store).

Let's create a shallow copy of the AntiKt4EMTopoJets and shift the pt of all jets up by 5%. In our source code, MyAnalysis/Root/MyxAODAnalysis.cxx, near the top add the include statement to do this shallow copy:

%SYNTAX{ syntax="cpp" }%

#include "xAODCore/ShallowCopy.h"

%ENDSYNTAX%

Now in execute(), we assume you already have some lines like:

%SYNTAX{ syntax="cpp" }%

xAOD::JetContainer::iterator jetSC_itr = (jets_shallowCopy.first)->begin();

xAOD::JetContainer::iterator jetSC_end = (jets_shallowCopy.first)->end();

for( ; jetSC_itr = jetSC_end; ++jetSC_itr ) {

// apply a shift in pt, up by 5%

double newPt = (*jetSC_itr)->pt() * (1 + 0.05);

// from: https://svnweb.cern.ch/trac/atlasoff/browser/Event/xAOD/xAODJet/trunk/xAODJet/versions/Jet_v1.h

xAOD::JetFourMom_t newp4(newPt, (*jetSC_itr)->eta(), (*jetSC_itr)->phi(), (*jetSC_itr)->m());

(*jetSC_itr)->setJetP4( newp4); // we've overwritten the 4-momentum

} // end iterator over jet shallow copy

%ENDSYNTAX%

By default when you create a shallowCopyContainer you take ownership of the pairs (meaning you need to take care to delete the both). You can give ownership to either the TStore or TEvent, which will then handle deletion for you. If you want to work with the pair internally to your algorithm or pass it around to different algorithms without writing it out to an output xAOD, you should give it to the TStore. If you plan to write it out to an output xAOD you can give ownership to TEvent. Each of these cases is described below:

( *jet_itr )->auxdecor< int >( "mySignal" ) = 1; // 1 = yes, it's a signal jet!

}

else{

( *jet_itr )->auxdecor< int >( "mySignal" ) = 0; // 0 = nope, not a signal jet

}

} // end for loop over jets

%ENDSYNTAX%

}

else{

( *jetSC_itr )->auxdata< int >( "mySignal" ) = 0; // 0 = nope, not a signal jet

}

} // end iterator over jet shallow copy

%ENDSYNTAX%

Here we have added an integer variable called 'mySignal' to the shallow copied jets

---+++ Future feature: Easy container slimming

As a final ingredient to writing out modified objects, you can select which of their properties should be written to the output file. The xAOD design was based around the idea that objects/containers may be slimmed during the analysis easily.

As you should know, all the data payload of xAOD objects/containers is in their auxiliary store objects. Because of this the way to specify which variables should be written out, is to set a property for the auxiliary store in question. This is done using the xAOD::TEvent::setAuxItemList function. By putting something like this into your algorithm's initialize() function:

%SYNTAX{ syntax="cpp" }%

// Set which variables not to write out:

event->setAuxItemList( "AntiKt4EMTopoJetsAux.", "-NumTrkPt1000.-NumTrkPt500" );

// Set which variable to do write out:

event->setAuxItemList( "GoodJetsAux.", "JetGhostArea.TrackCount" );

%ENDSYNTAX%

Unfortunately at the time of writing this still has some issues when using multiple input files (code crashing on a few files into the job), but the formalism is going to be this once the code works as intended...

---+ 10. More EventLoop features

!EventLoop is an ASG supported ROOT tool for handling optimized event looping and management. Below are some nice features of using EventLoop in your analysis. One nice feature is the different drivers available to easily switch between running your code locally, on the Grid, using PROOF, or on a batch system.

A more complete description of the tool can be found on the dedicated EventLoop wiki page.

---++ Running on the grid

In this section of the tutorial we will teach you how to run on the grid. There are two main advantages to running on the grid: First, you have access to the vast computing power of the grid, which means even very lengthy jobs can finish within hours, instead of days or weeks. Second, you have direct access to all the datasets available on the grid, which means you don't need to download them to your site first, which can save you hours if not days of time. There are also two main disadvantages: First, for all but the simplest jobs your turnaround time will be measured in hours, if not days. And this time can vary depending on the load at the various grid sites. Second, there are more things beyond your control that can go wrong, e.g. the only grid site with your samples may experience problems and go offline for a day or two thereby delaying the execution of your jobs.

If you need more information on options available for running on the grid, check out the grid driver documentation.

---++ Using TTreeCache

The speed with which your jobs will run will in most cases be dominated by the speed with which you can read your input data. So it is often worthwhile to try to maximize that speed.

One way to improve read performance is to use TTreeCache, which will predict the data that you are likely to access and preload it in big chunks. This is mostly important when reading files over the network, and it is virtually mandatory in case you read the files directly from the grid (see the next section on FAX).

Using TTreeCache with EventLoop is very straightforward, just specify the size of the cache you would like for your job before submitting it (in this case 10MB):

%SYNTAX{ syntax="cpp" }%

job.options()->setDouble (EL::Job::optCacheSize, 10*1024*1024);

%ENDSYNTAX%

The default way in which TTreeCache predicts what data your job will access is by looking at the first n events and assume that the pattern will hold for the rest of your jobs. If you want to, you can change the number of events to be read to suit your needs:

%SYNTAX{ syntax="cpp" }%

job.options()->setDouble (EL::Job::optCacheLearnEntries, 20);

%ENDSYNTAX%

You may have to play around a little to determine which number works best for you. There is a definite tradeoff, in that a too large number will mean that you read too many events without the benefit of the cache, while a too small number means that you will not be able to cache variables that you do not access on every event.

---++ Using FAX

FAX is a mechanism that allows you to read data directly from the grid without downloading it first. Depending on your usage pattern this may not only be more convenient, but also faster than downloading the files directly: If you download a file, you have to download it in its entirety, whereas if you read the file directly from the grid you may read only a small fraction of it. So if you are reading the file only once (or only a few times) this may be an option to improve your analysis workflow.

%W% Warning: While not strictly necessary, it is strongly recommended that you use TTreeCache when using FAX, otherwise your performance is likely to be very poor. So if you haven't done so already, you should work through the section above on TTreeCache first.

*Exit your shell and start a new session on lxplus (or your working environment).* As always, we have to setup the atlas environment:

%SYNTAX{ syntax="sh" }%

setupATLAS

%ENDSYNTAX%

Now in your new shell source the environments to use the FAX tools and establish a voms proxy:

%SYNTAX{ syntax="sh" }%

lsetup fax

voms-proxy-init -voms atlas

%ENDSYNTAX%

And when prompted enter your grid certificate password. The first line here takes care of setting up FAX and other necessary tools. If you do not work on lxplus you may need to do something different. The second line will initialize the VOMS proxy that holds your grid certificate. If you have problems with that command you should ask your resident grid expert. There are too many things that can go wrong with that to discuss them here.

Navigate to your working area, and from there setup your Analysis Release, following the recommendations above What to do everytime you log in.