DxAOD 解析

ATLAS-Japanソフトウェア講習会2016

東京大学@2016年12月27日

イントロダクション

このページは以下のATLASチュートリアルを元にしています。

https://twiki.cern.ch/twiki/bin/view/AtlasComputing/SoftwareTutorialxAODAnalysisInROOT?;


ATLAS実験では、AODと呼ばれるファイルに物理解析用のデータが保存されます。LHC-Run1ではAODはAthenaでしかアクセスすることが出来ませんでしたが、LHC-Run2ではROOTでも解析出来るようになりました。その際、AODのことをxAODと呼ぶようになりました。実際のデータ解析はこのxAODを用いて行うわけですが、膨大な衝突データの全てを見るのはサイズがデカすぎて不便です。そこで、 DerivationFramework でxAODに対して簡単な事象選別を行い、サイズを減らしたDxAODと呼ばれるファイルを用います。DxAODのフォーマットはxAODと全く同じです。xAODからDxAODに変換する際、簡単なバグ修正やエネルギー較正のアップデートなどを DerivationFramework を走らせるだけで行え、コンピュータリソースの無駄遣いを抑制することができます。したがって、Derivation前のxAODを直で読みにいくことは推奨されていません。どのバージョンの DerivationFramework で作られたDxAODかを常に確認する必要があります。そうした情報は以下のページにまとまっています。

https://twiki.cern.ch/twiki/bin/view/AtlasProtected/DerivationProductionTeam#Info_on_AtlasDerivation_caches_a

DxAOD が作られるまでに検出器で得られた信号は物理的に意味のある変数(運動量など)に変換されます。その際に種々の較正がかけられているのですが、全てを完璧にチューニングする事はできません。さらに詳細な較正を解析レベルで行います。また、MCシミュレーションが常に完璧でははないため、検出効率、運動量スケール、分解能といったパラメータの実データとのズレを補正する必要があります。Combined Perforamce(CP)グループがこれらの解析に対して推奨ツールを用意しており、CP toolと呼びます。 ROOTを用いる場合に、これらCP toolのチェックアウトやコンパイル、解析用のパッケージの生成などを行うビルドシステムが RootCore です。

https://twiki.cern.ch/twiki/bin/view/AtlasComputing/RootCore?redirectedfrom=Atlas.RootCore

このチュートリアルでは以下の内容を実際に手を動かしながら覚えていきます:

  • 新しいパッケージを作る方法
  • xAODへのアクセスの方法とCP toolの使い方
  • xAODからHistogramやNTupleを作る方法

環境設定

作業はlogin.icepp.jpで行います。

まず、作業用のディレクトリを作成しましょう。

$ mkdir ROOTAnalysisTutorial
$ cd ROOTAnalysisTutorial

以下のコマンドを実行してAnalysis Releaseが使えるようにします。Analysis Releaseを設定することで、 RootCore, CP tool等が使えるようになります。

$ setupATLAS 
$ rcSetup Base,2.4.23
$ rc find_packages
$ rc compile

今回は2.4.23というバージョンを設定しました。Analysis Releaseのバージョンによって、 CP toolのバーションも管理されています。なので、Analysis Releaseのバージョンさえ確認しておけば、 解析グループによってうっかり異なるバージョンのCP toolを使ってしまうことがないようになっています。 以降、新しいセッションを始める場合には以下のコマンドのみで設定ができます。

$ setupATLAS
$ rcSetup

各releaseでどのパッケージにどのような変更があるかは

https://twiki.cern.ch/twiki/bin/viewauth/AtlasProtected/AnalysisBaseReleaseNotes24

にまとまっています。また、

$ rc version

で各パッケージのversionを知ることができます。

Releaseを変更する場合は

$ rcSetup -u
$ rcSetup Base,X.X.XX
$ rc find_packages
$ rc compile

です。

データ

せっかくなのでATLASの実データを触っていきましょう。

最新の1run分のDxAODを以下のdirectoryにコピーしておきました。

/gpfs/fs2001/nobe/data2086b/sample_forTutorial2016/data16_13TeV.00311481.physics_Main.merge.DAOD_HIGG2D4.f758_m1710_p2880/

今回はHiggs解析グループで主に使用されるHIGG2D4という種類のDxAODを用います。DxAODを作る際に予めかけられている事象選別は以下に簡単にまとまっています。

https://docs.google.com/spreadsheets/d/1vzqMXKpu5yuXGAX0Pbb6NmY35aQa_jU_sj0kFU29mwI/edit#gid=1484460466

2本のレプトンが要求されているため、HIGG2D4を用いてレプトンを要求しない解析(例えば、tt -> bqqbqqなど)を行うことはできません。自分の目的の物理に応じた種類のDxAODを使うようにしましょう。

自分の解析用パッケージの作成

それでは、新しくパッケージを作って解析をしてみます。

新規パッケージの作成

自分のに作業ディレクトリに移動してパッケージを作ります。

$ rc make_skeleton MyNewPackage

これでほぼ空っぽの MyNewPackage が作られました。

$ ls MyNewPackage
MyNewPackage Root cmt

アルゴリズムの作成

RootCore で提供されている EvenrLoop パッケージを使ってアルゴリズムクラスを作成します。

$ $ROOTCOREBIN/user_scripts/EventLoop/make_skeleton MyNewPackage MyAnalysisAlg

以下のファイルが作成されました。initialize(),finalize(),execute()などが定義されており、用途が推測できます。

この MyNewPackage をコンパイルします。

$ rc find_packages
$ rc compile

MyNewPackage が認識されてコンパイルされました。それではこれらのファイルを編集していきましょう。 xAODRootAccess というパッケージを介して、xAODにアクセスします。 MyAnalysisAlg.cxxを編集していきます。まずヘッダーファイルを追加します。


// Infrastructure include(s): 
#include "xAODRootAccess/Init.h"
#include "xAODRootAccess/TEvent.h" // ASG status code check
#include <AsgTools/MessageCheck.h>

次に、setupJob(EL::Job& job)に以下を追加します:

job.useXAOD (); ANA_CHECK(xAOD::Init());<br />

次にinitialize()に以下を追加します。

xAOD::TEvent* event = wk()-&gt;xaodEvent(); Info("initialize()", "Number of events = %lli", event-&gt;getEntries() );<br />

MyNewPackage /cmt/Makefile.RootCoreに依存するパッケージをリストします:

PACKAGE_DEP = EventLoop xAODRootAccess AsgTools <br />

コンパイルします。

$ rc find_packages
$ rc compile

解析コードの本体を作成

解析コードの本体を作成します。

$ mkdir MyNewPackage/util
$ touch MyNewPackages/util/testRun.cxx

testRun.cxxは以下の内容にします。


<br />#include "xAODRootAccess/Init.h" <br />#include "SampleHandler/SampleHandler.h" <br />#include "SampleHandler/ScanDir.h" <br />#include "SampleHandler/ToolsDiscovery.h" <br />#include "EventLoop/Job.h" <br />#include "EventLoop/DirectDriver.h" <br />#include "SampleHandler/DiskListLocal.h" <br />#include

<br />#include "MyNewPackage/MyAnalysisAlg.h"

int main( int argc, char* argv[] ) {

 // Take the submit directory from the input if provided: std::string submitDir = "submitDir"; if( argc &gt; 1 ) submitDir = argv[ 1 ];

 // Set up the job for xAOD access: xAOD::Init().ignore();

 // Construct the samples to run on: SH::SampleHandler sh;

 // use SampleHandler to scan all of the subdirectories of a directory for particular MC single file: const char* inputFilePath = gSystem-&gt;ExpandPathName ("/gpfs/fs2001/nobe/data2086b/sample_forTutorial2016/data16_13TeV.00311481.physics_Main.merge.DAOD_HIGG2D4.f758_m1710_p2880/"); SH::ScanDir().filePattern("*pool.root*").scan(sh,inputFilePath);

 // Set the name of the input TTree. It's always "CollectionTree" // for xAOD files. sh.setMetaString( "nc_tree", "CollectionTree" );

 // Print what we found: sh.print();

 // Create an EventLoop job: EL::Job job; job.sampleHandler( sh ); job.options()-&gt;setDouble (EL::Job::optMaxEvents, 500); // maximum number of events to be analyzed

 // Add our analysis to the job: MyAnalysisAlg* alg = new MyAnalysisAlg (); job.algsAdd( alg );

 // Run the job using the local/direct driver: EL::DirectDriver driver; driver.submit( job, submitDir );

 return 0; }<br />

コンパイルします。

$ rc find_packages
$ rc compile

以下のコマンドでtestRunというプログラムを走らせることができます。

$ testRun submitDir

submitDirはアウトプットファイルが出来る場所で、既にあるディレクトリ名を指定するとジョブは走りません。

このままだとEventLoopで何もしないで最初の500eventsを読み込むだけのコードです。

データにアクセスする方法

以下のコマンドでAODにどのようなデータが保存されているか確認出来ます。

$ checkxAOD.py /gpfs/fs2001/nobe/data2086b/sample_forTutorial2016/data16_13TeV.00311481.physics_Main.merge.DAOD_HIGG2D4.f758_m1710_p2880/DAOD_HIGG2D4.09977266._000001.pool.root.1 

名前からどういうデータが保存されているか推測できます。 さらに、例えばMuons(DataVector<xaod::muon_v1>)というコンテナーにはどういう変数が保存されているかを以下のコマンド見ることができます。

?$ root -l /gpfs/fs2001/nobe/data2086b/sample_forTutorial2016/data16_13TeV.00311481.physics_Main.merge.DAOD_HIGG2D4.f758_m1710_p2880/DAOD_HIGG2D4.09977266._000001.pool.root.1 
root [1] CollectionTree->Print("Muons*")
pT, eta, phiなどに加えて様々な変数が定義されています。

データにアクセス (1)

Run NumberやEvent Numberにアクセスしてみます。MyAnalysisAlg.cxxを編集していきます。まず、ヘッダーファイルをインクルードします:

%CODE{"cpp"}%

// EDM includes:

#include "xAODEventInfo/EventInfo.h

$ENDCODE%

execute()に以下を追加します:

xAOD::TEvent* event = wk()-&gt;xaodEvent();<br /> <br />// Event information<br />const xAOD::EventInfo* eventInfo = 0;<br />ANA_CHECK(event-&gt;retrieve(eventInfo, "EventInfo")); <br /> <br />int runNumber  = eventInfo-&gt;runNumber();<br />int eventNumber = eventInfo-&gt;eventNumber();<br /> <br />Info("execute()", "Run = %i : Event = %i", runNumber, eventNumber);

RootCoreにxAODEventInfoを加えます:

%CODE{"cpp"}%

PACKAGE_DEP = EventLoop xAODRootAccess AsgTools xAODEventInfo

%ENDCODE%

コンパイルした後、プログラムを実行してみます。 $ rc find_packages $ rc compile $ rm -r submitDir/ $ testRun submitDir ... Running sample: mc15_13TeV MyAnalysisAlg::initial... INFO Number of events = 75718 MyAnalysisAlg::execute() INFO Run = 222526 : Event = 1887985 MyAnalysisAlg::execute() INFO Run = 222526 : Event = 1887757 MyAnalysisAlg::execute() INFO Run = 222526 : Event = 1886250 ...

Run NumberとEvent Numberにアクセスすることが出来ました。

さて、 ここでGoodRunsListSelectionToolというCP toolを導入してみましょう。 GoodRunsListSelectionToolはxAODEventInfoを入力にします。 GoodRunListというのは物理解析にとって良いルミブロックであるかどうかがxmlで記述されたものです。 つまりATLASで検出器の状態が良くなかったと判断されたルミブロックをこのツールで除いていくことが出来ます。


-- TakuyaNobe - 2016-12-06

Edit | Attach | Watch | Print version | History: r5 < r4 < r3 < r2 < r1 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r5 - 2016-12-09 - TakuyaNobe
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    Main All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2019 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback