TWiki> Main Web>TWikiUsers>MasahiroMorinaga>TokyoICEPPATLASTutorialPage>DataAndMc (2014-04-24, MasahiroMorinaga) EditAttachPDF
-- MasahiroMorinaga - 11 Apr 2014

ここでは、実際にデータやモンテカルロを使って解析を行うときに必要なことを示します。ATLAS analysisのページだけでは不十分です。

シミュレーションが完璧ではないのでモンテカルロとデータは一致しない場合がほとんどです、例えばトリガーエフィシェンシー等。それらを合わせなければデータとバックグランドはきれいに合いません。

そのため、各パフォーマンスグループがせっせといろんなツールを作ってくれています、ここではこれを使えるようになりましょう。(SUSYグループにはSUSYToolというゆとりつーるがあるみたいですが、一度は自分で正しく使えるようになりましょう。)

ATLAS Official Recommendations

以下に、ATLASのオフィシャルなrecommendationのツールをおいておきますので、svnからとってきてください。

使い方は個々のMakefileやページにかいてあります、使い方がわからなければ遠慮せずに聞いてください。

SVN packages:

[[https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/D3PDTools/RootCore/tags/RootCore-00-01-95] [https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/D3PDTools/RootCore/tags/RootCore-00-01-95]]
https://svnweb.cern.ch/trac/atlasoff/browser/DataQuality/GoodRunsLists/tags/GoodRunsLists-00-01-06
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/AnalysisCommon/PileupReweighting/tags/PileupReweighting-00-02-12
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/TileID/TileTripReader/tags/TileTripReader-00-00-19
https://svnweb.cern.ch/trac/atlasoff/browser/Trigger/TrigAnalysis/TrigRootAnalysis/tags/TrigRootAnalysis-00-01-05
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/egamma/egammaAnalysis/egammaAnalysisUtils/tags/egammaAnalysisUtils-00-04-55
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/tags/ElectronEfficiencyCorrection-00-00-50
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/AnalysisCommon/PATCore/tags/PATCore-00-00-16
https://svnweb.cern.ch/trac/atlasoff/browser/Trigger/TrigAnalysis/TrigMuonEfficiency/tags/TrigMuonEfficiency-00-02-48
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonMomentumCorrections/tags/MuonMomentumCorrections-00-09-23
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonEfficiencyCorrections/tags/MuonEfficiencyCorrections-02-01-19
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonIsolationCorrection/tags/MuonIsolationCorrection-01-01
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/JetTagging/JetTagPerformanceCalibration/CalibrationDataInterface/tags/CalibrationDataInterface-00-03-06
https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/JetMissingEtID/JetSelectorTools/tags/JetSelectorTools-00-01-01
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/Jet/ApplyJetCalibration/tags/ApplyJetCalibration-00-03-20
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/Jet/JetUncertainties/tags/JetUncertainties-00-08-07
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/Jet/JetResolution/tags/JetResolution-02-00-02
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/Jet/JetAnalysisTools/ApplyJetResolutionSmearing/tags/ApplyJetResolutionSmearing-00-01-02
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/Jet/JetAnalysisTools/JVFUncertaintyTool/tags/JVFUncertaintyTool-00-00-03
https://svnweb.cern.ch/trac/atlasoff/browser/Reconstruction/MissingETUtility/tags/MissingETUtility-01-02-05

Twiki pages:

前段階

解析をする前段階としてしなければいけないことがいくつかあります、大きく分けると

データの取得状態が良いものを選ぶ またはモンテカルロを補正する

データ取得がうまく言っているイベントを選びましょう。判断基準は複数あります。

・ GRL(Good Runs List) : 検出器のパフォーマンスが悪い状況で取得したデータは解析には使えません。これはGRL ToolというものとGRL XMLファイルというもので選んでくれます。

・ Pileup reweighting : モンテカルロに入っているパイルアップコンディションは実際とは違います、これを補正します。(スケールファクターをかける)

・ Event Cleaning : GRL等では対処できない、個々の悪い状況でのイベントを消してやります。例えば、カロリメータでの測定がうまくいっていないイベント等。 


//GRL
if(m_GoodRunsList.HasRunLumiBlock(TheRunNumber,lbn))

// PU reweight
m_tPileUpTool->SetRandomSeed(mc_channel_number+EventNumber);
int TheRunNumber      = m_tPileUpTool->GetRandomRunNumber(RunNumber);
int AverageIntPerXing = (lbn==1 && int(averageIntPerXing+0.5)==1) ? 0. : averageIntPerXing;
int WeightCombined    = m_tPileUpTool->GetCombinedWeight(RunNumber,mc_channel_number,AverageIntPerXing);

// Cleaning
if(larError < 2 && tileError!=2 && ((coreFlags&0x40000) == 0 && m_TileReader->checkEvent(RunNumber,lbn,EventNumber)!=0))

見たい物理事象を効率よく選ぶ

解析をする上で効率よくイベントを選みましょう。

トリガー

自分の見たい物理事象に合ったトリガーが鳴ったイベントを見ましょう。このWW->lvqqの場合、レプトンはひとつだけなので、シングルレプトントリガーを使います。single-lepton triggerはATLASで一番優先されるトリガーです。簡単にまとめると

チャンネル トリガー しきい値 その他
Electron EF_e24vhi_medium1 pT>24GeV mediumなIDとアイソレーションがかかっています
Electron EF_e60_medium1 pT>60GeV mediumなIDがかかっています
Muon EF_mu24i_tight pT>24GeV tightなIDとアイソレーションがかかっています
Muon EF_mu36_tight pT>36GeV mediumがかかっています

となっています、それぞれElectron/Muonチャンネルで2つのトリガーのORを取りましょう。

・トリガーマッチング : オンラインでトリガーした物理オブジェクトがオフラインで再構成されたオブジェクトとマッチしているかをチェックします。TrigMatching Toolを使いましょう。 ・トリガーエフィシェンシー : トリガーエフィシェンシーはデータとモンテカルロで合っていないことがありますので、それをスケールファクターを使って補正します。Electron scale factor toolやMuon Trigger scale factor toolを使いましょう。

Primary Vertex(PV)とか

PVにアソシエイトしているトラックが3本以上あることを要求します。 

if(vxp_nTracks->at(0) >= 3)

各補正やセレクションの再計算

Electronの補正

使うツール 使い方 data or mc 説明
方向 特になし ElectronsEnergyDirection どっちも トラッカーのヒットが多いものはトラックの方向、そうじゃないものはクラスターのを使う
エネルギー egammaAnalysisUtils 下の dataもmcも ntuplizeした時点からアップデートとかされるので。。
ID egammaAnalysisUtils 下だよ dataもmcも バグがあったり改善するので
Isolation CaloIsolationCorrections 下ね どっちも 普通の解析はパイルアップの補正を追加でします。
IDのSF TElectronEfficiencyCorrectionToolとPATCore 下ですね mcだけ mcは完璧じゃないのでefficiencyの違いを補正します。

下にコードの例を示しておきます、ざっくりしているのでコンパイルする際は、じぶんで！適宜変更してください。 

//direction
int nHit = (el_nPixHits->at(i)+el_nSCTHits->at(i))
float eta,phi;
if(nHit < 4){ eta = el_cl_eta->at(i); phi = el_cl_phi->at(i); } else { eta = el_tracketa->at(i); phi = el_trackphi->at(o); }

//Energy 
float e;
if(isMC){
  m_egammaTool->SetRandomSeed(seed);
  double smearcorr = m_egammaTool->getSmearingCorrection(cl_eta,cl_E,egRescaler::EnergyRescalerUpgrade::NOMINAL);
  e = cl_E * smearcorr;
}else{
  e = m_egammaTool->applyEnergyCorrection(cl_eta,cl_E,egRescaler::EnergyRescalerUpgrade::Electron,egRescaler::EnergyRescalerUpgrade::Nominal);
}

// ID
bool LoosePP = isLoosePlusPlus(eta,eT,rHad,rHad1,Reta,w2,f1,wstot,DEmaxs1,deltaEta,nSi,nSiOutliers,nPix,nPixOutliers,egammaMenu::eg2012,false,false);
bool MediumPP = isMediumPlusPlus(eta,eT,f3,rHad,rHad1,Reta,w2,f1,wstot,DEmaxs1,deltaEta,d0,TRratio,nTRT,nTRTOutliers,nSi,nSiOutliers,nPix,nPixOutliers,nBlayer,nBlayerOutliers,expectBlayer,egammaMenu::eg2012,false,false);
bool TightPP = isTightPlusPlus(eta,eT,f3,rHad,rHad1,Reta,w2,f1,wstot,DEmaxs1,deltaEta,d0,TRratio,nTRT,nTRTOutliers,nSi,nSiOutliers,nPix,nPixOutliers,nBlayer,nBlayerOutliers,eT,eOverp,deltaphi,convBit,egammaMenu::eg2012,false,false);

// Isolation
float Et20 = CaloIsoCorrection::GetPtNPVCorrectedIsolation(nPV2,energy,etaS2,etaPointing,etaCluster,20,isMC,Etcone_value20,false,CaloIsoCorrection::ELECTRON);

// ID SF
PATCore::ParticleDataType::DataType data_type = (isAF2)?  PATCore::ParticleDataType::Fast : PATCore::ParticleDataType::Full;
// reco SF
const Root::TResult &result_reco = m_el_reco_SF->calculate(data_type,runNumber,eta,pt);
// tightPP
const Root::TResult &result_ID = m_el_tight_SF->calculate(data_type,runNumber,eta,pt);
float sf = result_ID.getScaleFactor() * result_reco.getScaleFactor();

Muonの補正

Muonは結構簡単です。

使うツール 使い方 data or mc 説明
エネルギー MuonMomentumCorrections 下の mcだけ データにあわせてやります
Isolation MuonIsolationCorrection 下ね どっちも?? パイルアップの補正を追加でします。
IDのSF AnalysisMuonConfigurableScaleFactors 下ですね mcだけ mcは完璧じゃないのでefficiencyの違いを補正します。 

// Energy
float ptcb = mu_staco_pt;
float ptid = (mu_staco_id_qoverp_exPV->at(i) != 0.) ? fabs(sin(mu_staco_id_theta_exPV->at(i))/mu_staco_id_qoverp_exPV->at(i)) : 0.;
float ptms = (mu_staco_me_qoverp_exPV->at(i) != 0.) ? fabs(sin(mu_staco_me_theta_exPV->at(i))/mu_staco_me_qoverp_exPV->at(i)) : 0.;
float corr_pt = ptcb;
if(isMC) {
  m_muonTool->SetSeed(EventNumber, i);
  m_muonTool->Event(ptms,ptid,ptcb,eta,charge);
  corr_pt = m_muonTool->pTCB();
}

// Isolation
float et20_corr = m_muisoTool->CorrectEtCone(etcone20,nPV_3track,eta,"cone20Comb"))


// ID SF
float sf = m_stacoSF->scaleFactor(charge,tlv);

Jetの補正

ジェットの補正はかいせきに結構効きます。 まあ基本的には4vectorのキャリブレーションです。 

    Jet_E         = jet_AntiKt4LCTopo_constscale_E->at(i);
    Jet_M         = jet_AntiKt4LCTopo_constscale_m->at(i);
    Jet_Eta       = jet_AntiKt4LCTopo_constscale_eta->at(i);
    Jet_Phi       = jet_AntiKt4LCTopo_constscale_phi->at(i);
    Jet_Ax        = jet_AntiKt4LCTopo_ActiveAreaPx->at(i);
    Jet_Ay        = jet_AntiKt4LCTopo_ActiveAreaPy->at(i);
    Jet_Az        = jet_AntiKt4LCTopo_ActiveAreaPz->at(i);
    Jet_Ae        = jet_AntiKt4LCTopo_ActiveAreaE->at(i);
    Jet_rho       = Eventshape_rhoKt4LC;

    TLorentzVector ajet = tool->m_jetCalibTool->ApplyJetAreaOffsetEtaJES(Jet_E,Jet_Eta,Jet_Phi,Jet_M,Jet_Ax,Jet_Ay,Jet_Az,Jet_Ae,Jet_rho,EventInfo->Mu(),nPV_2trks);

あ、あとb-tagのSFも計算しましょう。mcだけです。ややこしいので少しだけ書いておきます、細かいところは自分でCalibrationDataInterfaceを見て調べてね。 

  float sf = 1.0;
  bool isJVFOK = (fabs(EMJet_Eta)<2.4 && fabs(JVF)>0.5 && EMJet_Pt<50.);
  if(fabs(EMJet_Eta) < 2.5){
    if(isJVFOK){ BtagJet->jetAuthor = "AntiKt4TopoLCJVF0_5"; }
    else{        BtagJet->jetAuthor = "AntiKt4TopoLCnoJVF";  }

    if(fabs(EMJet_Eta) < 2.5 && EMJet_Pt > 20.){
      std::string flavor = "";
      if     (flavor_truth_label == 4){ flavor = "C"; }
      else if(flavor_truth_label == 5){ flavor = "B"; }
      else if(flavor_truth_label ==15){ flavor = "T"; }
      else                            { flavor = "Light"; }

      if(MV1>0.7892){
        Analysis::CalibResult TagResult = BtagCalib->getScaleFactor((*BtagJet),flavor,"0_7892",BtagUncertainty);
        sf = TagResult.first;
      }else{
        Analysis::CalibResult NotTagResult = BtagCalib->getInefficiencyScaleFactor((*BtagJet),flavor,"0_7892",BtagUncertainty);
        sf = NotTagResult.first;
      }

METの補正

METの計算は実はややこしいというのはすでに言いましたが、以下に詳細を示しておきます。基本的には、各物理オブジェクトのエネルギーや方向の補正が行われたら、その効果をMETの計算に伝搬させてやります。

これにはMissingETUtilityというパッケージを用います。 

  vector<vector<float> > new_wet, new_wpx, new_wpy; new_wet.clear(); new_wpx.clear(); new_wpy.clear();
  vector<vector<unsigned int> > new_statusWord; new_statusWord.clear();
  for(int i=0;i<mu_staco_n;i++) {
    vector<float> univec;
    univec.push_back(1.);
    new_wet.push_back(univec);
    new_wpx.push_back(univec);
    new_wpy.push_back(univec);
    vector<unsigned int> defvec;
    defvec.push_back(MissingETTags::DEFAULT);
    new_statusWord.push_back(defvec);
  }

  tool->m_METUtil->reset();
  tool->m_METUtil->setJetParameters(jet_smeared_pt,jet_AntiKt4LCTopo_eta,jet_AntiKt4LCTopo_phi,jet_AntiKt4LCTopo_E,jet_AntiKt4LCTopo_MET_wet,jet_AntiKt4LCTopo_MET_wpx,jet_AntiKt4LCTopo_MET_wpy,jet_AntiKt4LCTopo_MET_statusWord);
  tool->m_METUtil->setOriJetParameters(jet_AntiKt4LCTopo_pt);
  tool->m_METUtil->setElectronParameters(el_smeared_pt,el_smeared_eta,el_smeared_phi,el_MET_wet,el_MET_wpx,el_MET_wpy,el_MET_statusWord);
  tool->m_METUtil->setMuonParameters(mu_smeared_pt,mu_staco_eta,mu_staco_phi,&new_wet,&new_wpx,&new_wpy,&new_statusWord);
  tool->m_METUtil->setExtraMuonParameters(mu_staco_ms_qoverp,mu_staco_ms_theta,mu_staco_ms_phi,mu_staco_charge);
  tool->m_METUtil->setMETTerm(METUtil::RefGamma,MET_RefGamma_etx_CentralReg,MET_RefGamma_ety_CentralReg,MET_RefGamma_sumet);
  tool->m_METUtil->setMETTerm(METUtil::RefTau,MET_RefTau_etx_CentralReg,MET_RefTau_ety_CentralReg,MET_RefTau_sumet);
  tool->m_METUtil->setMETTerm(METUtil::CellOutEflow,MET_CellOut_Eflow_et * cos(MET_CellOut_Eflow_phi),MET_CellOut_Eflow_et * sin(MET_CellOut_Eflow_phi),MET_CellOut_Eflow_sumet);

  
  METObject RefFinal     = tool->m_METUtil->getMissingET(METUtil::RefFinal);
  METObject RefEle       = tool->m_METUtil->getMissingET(METUtil::RefEle);
  METObject RefGamma     = tool->m_METUtil->getMissingET(METUtil::RefGamma);
  METObject RefJet       = tool->m_METUtil->getMissingET(METUtil::RefJet);
  METObject RefTau       = tool->m_METUtil->getMissingET(METUtil::RefTau);
  METObject RefMuon      = tool->m_METUtil->getMissingET(METUtil::RefMuon);
  METObject CellOutEflow = tool->m_METUtil->getMissingET(METUtil::CellOutEflow);
  METObject MuonTotal    = tool->m_METUtil->getMissingET(METUtil::MuonTotal);

el/mu/jetに関しては、smearedというvectorを使っています、これは上で補正やキャリブレーションしなおした際の値です。すべてのel/mu/jetに関してどこかに保存しておいてください。

ここまでで、どうしてもわからないこと等があるひとは、森永のコード覗いてもらっても良いです（どっちにしても、よくわからないと思いますが。）/home/morinaga/vvReso/2014Mar/Code 以下にあります。

きちんと補正できているかのチェック

各種補正をしてきましたが、自分の書いたコードが正しく動いているかのチェックはとっっっっっっっても大切です。解析の結果を出すときは他の人と独立したコードでAcceptance Challengeということをします。ACとか呼びます。 具体的には同じデータセットを使って、el/mu/jet等のオブジェクトのセレクションをすべての段階で合計数をあわせます。もちろん本来のカットフローもスケールファクター込みで合わせます。

この作業は結構つらい作業です。(SUSYToolsなどのゆとりツールを使っていれば楽ですが。。。)

というわけで、ここではそれを行います。下に森永のカットフローを載せておくので自分のものとチェックしてください。 都合上、 

/gpfs/fs3001/atljphys/datafiles/group.phys-exotics/mc12_8TeV/back/user.matkinso.mc12_8TeV.158234.CalcHepPythia8_AU2CTEQ6L1_KKGravitonWW_lvjj_m1200_1LeptonPt8.LVQQ2.e1370_s1499_s1504_r3658_r3549/*.root
というサンプルを使いますが、動かない人は自分のコード直してね。(あまりまじめにやってないので森永が間違っている可能性をあるので気軽に聞いてください。)

LeadingというのはPtが最も大きいものという意味です。step3のコードに 

sort(GoodJets.begin(), GoodJets.end(), [](const TLorentzVector l,const TLorentzVector m) -> bool { return l.Pt() > m.Pt(); } );
というのを追加しておいたので参考にしてください。（Ptの大きい順にソートしただけです。）

Selection 数
Electron
All 180244.0
Author 79902.0
Eta 77644.0
Charge 77644.0
ID 6011.0
OQ 5996.0
Pt 5621.0
Z0 5619.0
D0 5598.0
TrkIso 5589.0
CaloIso 5570.0
Muon
All 13578.0
Author 9552.0
Eta 9333.0
Combined 9333.0
Blayer 9127.0
Pix 8961.0
SCT 8942.0
Holes 8942.0
TRT 8934.0
Pt 6704.0
Z0 6697.0
D0 6521.0
QoverP 6438.0
TrkIso 6259.0
CaloIso 6189.0
AntiKt4 Jet
All 140568.0
Pt 37557.0
Eta 35370.0
JVF 34029.0
Cleaning 34025.0

個別のWeightもチェックしましょう。

EventNumber Channel(Electron: 1, Muon: 2) TrigSF PU Weight Lepton Pt(GeV) Lepton SF MET(GeV) Leading AntiKt4 Jet Pt(GeV)
10001 2 0.9639 2.1868 27.0909 0.9973 653.2912 734.6100
10002 1 0.9871 2.1868 182.0027 0.9666 293.8213 487.0157
10004 1 0.9872 2.1868 200.0938 0.9966 378.6221 555.8614
10007 2 0.9962 2.1868 515.0761 0.9973 73.5905 544.0573
10008 2 0.9754 2.1868 264.4872 0.9980 315.5975 597.2946
10009 2 0.9970 2.1868 175.9823 0.9932 277.2847 440.8742
10010 1 0.9837 2.1868 227.5870 1.0312 243.8748 466.1188
10011 1 0.9901 2.1868 202.2779 0.9666 258.8592 418.3131
10012 1 0.9926 2.1868 130.1272 0.9936 417.1374 549.7065
10013 1 1.0055 2.1868 38.9232 0.9703 524.7440 597.7899
10014 1 0.9898 2.1868 238.4951 0.9936 335.2587 666.3776
10015 1 0.9932 2.1868 293.3914 0.9666 220.2968 573.2068
10018 1 0.9890 2.1868 252.1191 0.9833 342.9873 584.1941
10019 2 0.9274 2.1868 349.5846 0.9941 232.8985 584.8477
10020 1 0.9932 2.1868 364.1934 0.9666 274.1742 594.2462
10021 2 1.0890 2.1868 53.0992 0.9997 221.8162 266.4462
10022 1 0.9815 2.1868 151.7365 0.9657 471.6005 622.1624
10024 2 0.9502 2.1868 190.9913 0.9984 252.0244 539.5239
10025 1 0.9873 2.1868 251.9504 0.9936 260.6291 507.4990
10026 2 1.0061 2.1868 160.4221 1.0011 482.4016 589.0309
10027 1 0.9932 2.1868 281.1481 0.9666 199.2651 397.8934

だいたい、あったら次のページにうつってください。

Edit | Attach | Watch | Print version | History: r6 < r5 < r4 < r3 < r2 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r6 - 2014-04-24 - MasahiroMorinaga
 
    • 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-2024 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
or Ideas, requests, problems regarding TWiki? use Discourse or Send feedback