Update and addition of LS analysis (Renu, Giacomo, Francesco)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisTaskSEDplus.cxx
CommitLineData
d486095a 1/**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/////////////////////////////////////////////////////////////
17//
18// AliAnalysisTaskSE for the extraction of signal(e.g D+) of heavy flavor
19// decay candidates with the MC truth.
4afc48a2 20// Authors: Renu Bala, bala@to.infn.it
21// F. Prino, prino@to.infn.it
22// G. Ortona, ortona@to.infn.it
d486095a 23/////////////////////////////////////////////////////////////
24
25#include <TClonesArray.h>
26#include <TNtuple.h>
27#include <TList.h>
1f4e9722 28#include <TString.h>
d486095a 29#include <TH1F.h>
4afc48a2 30#include <TDatabasePDG.h>
b557eb43 31
32#include "AliAnalysisManager.h"
33#include "AliAODHandler.h"
d486095a 34#include "AliAODEvent.h"
35#include "AliAODVertex.h"
36#include "AliAODTrack.h"
37#include "AliAODMCHeader.h"
38#include "AliAODMCParticle.h"
39#include "AliAODRecoDecayHF3Prong.h"
40#include "AliAnalysisVertexingHF.h"
41#include "AliAnalysisTaskSE.h"
42#include "AliAnalysisTaskSEDplus.h"
43
44ClassImp(AliAnalysisTaskSEDplus)
45
46
47//________________________________________________________________________
48AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus():
49AliAnalysisTaskSE(),
50fOutput(0),
4afc48a2 51fHistNEvents(0),
d486095a 52fNtupleDplus(0),
4afc48a2 53fHistLS(0),
54fHistLStrip(0),
55fHistLStripcut(0),
56fHistOS(0),
57fHistOSbkg(0),
58fHistLSnoweight(0),
59fHistLSsinglecut(0),
1f4e9722 60fFillNtuple(kFALSE),
4afc48a2 61fReadMC(kFALSE),
62fDoLS(kFALSE),
d486095a 63fVHF(0)
64{
65 // Default constructor
d486095a 66}
67
68//________________________________________________________________________
1f4e9722 69AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus(const char *name,Bool_t fillNtuple):
d486095a 70AliAnalysisTaskSE(name),
4afc48a2 71fOutput(0),
72fHistNEvents(0),
d486095a 73fNtupleDplus(0),
4afc48a2 74fHistLS(0),
75fHistLStrip(0),
76fHistLStripcut(0),
77fHistOS(0),
78fHistOSbkg(0),
79fHistLSnoweight(0),
80fHistLSsinglecut(0),
1f4e9722 81fFillNtuple(fillNtuple),
4afc48a2 82fReadMC(kFALSE),
83fDoLS(kFALSE),
d486095a 84fVHF(0)
85{
86 // Default constructor
87
88 // Output slot #1 writes into a TList container
89 DefineOutput(1,TList::Class()); //My private output
1f4e9722 90
91 if(fFillNtuple){
92 // Output slot #2 writes into a TNtuple container
93 DefineOutput(2,TNtuple::Class()); //My private output
1f4e9722 94 }
d486095a 95}
96
97//________________________________________________________________________
98AliAnalysisTaskSEDplus::~AliAnalysisTaskSEDplus()
99{
100 // Destructor
101 if (fOutput) {
102 delete fOutput;
103 fOutput = 0;
104 }
4afc48a2 105 if(fNtupleDplus){
106 delete fNtupleDplus;
107 fNtupleDplus=0;
108 }
d486095a 109 if (fVHF) {
110 delete fVHF;
111 fVHF = 0;
112 }
4afc48a2 113
d486095a 114}
115
116//________________________________________________________________________
117void AliAnalysisTaskSEDplus::Init()
118{
119 // Initialization
120
121 if(fDebug > 1) printf("AnalysisTaskSEDplus::Init() \n");
122
123 gROOT->LoadMacro("ConfigVertexingHF.C");
124
125 fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()");
126 fVHF->PrintStatus();
127
128 return;
129}
130
131//________________________________________________________________________
132void AliAnalysisTaskSEDplus::UserCreateOutputObjects()
133{
134 // Create the output container
135 //
136 if(fDebug > 1) printf("AnalysisTaskSEDplus::UserCreateOutputObjects() \n");
137
138 // Several histograms are more conveniently managed in a TList
139 fOutput = new TList();
140 fOutput->SetOwner();
1f4e9722 141 fOutput->SetName("OutputHistos");
142
143 Int_t nPtBins=4;
144 TString hisname;
145 for(Int_t i=0;i<nPtBins;i++){
146 hisname.Form("hMassPt%d",i);
147 TH1F* hm=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
148 hm->Sumw2();
149 hisname.Form("hSigPt%d",i);
150 TH1F* hs=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
151 hs->Sumw2();
152 hisname.Form("hBkgPt%d",i);
153 TH1F* hb=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
154 hb->Sumw2();
155
156 hisname.Form("hMassPt%dTC",i);
157 TH1F* hmtc=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
158 hmtc->Sumw2();
159 hisname.Form("hSigPt%dTC",i);
160 TH1F* hstc=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
161 hstc->Sumw2();
162 hisname.Form("hBkgPt%dTC",i);
163 TH1F* hbtc=new TH1F(hisname.Data(),hisname.Data(),100,1.765,1.965);
164 hbtc->Sumw2();
165
166
167 fOutput->Add(hm);
168 fOutput->Add(hs);
169 fOutput->Add(hb);
170 fOutput->Add(hmtc);
171 fOutput->Add(hstc);
172 fOutput->Add(hbtc);
173 }
4afc48a2 174 fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5);
175 fHistNEvents->Sumw2();
176 fHistNEvents->SetMinimum(0);
177 fOutput->Add(fHistNEvents);
178 if(fDoLS){
179 Double_t massDplus=TDatabasePDG::Instance()->GetParticle(411)->Mass();
8c34bd86 180
4afc48a2 181 fHistLS = new TH1F("fHistLS","fHistLS",100,massDplus-0.2,massDplus+0.2);
182 fHistLS->Sumw2();
183 fOutput->Add(fHistLS);
184
185 fHistLStrip = new TH1F("fHistLStrip","fHistLStrip",100,massDplus-0.2,massDplus+0.2);
186 fHistLStrip->Sumw2();
187 fOutput->Add(fHistLStrip);
188
189 fHistLStripcut = new TH1F("fHistLStripcut","fHistLStripcut",100,massDplus-0.2,massDplus+0.2);
190 fHistLStripcut->Sumw2();
191 fOutput->Add(fHistLStripcut);
192
193 fHistOS = new TH1F("fHistOS","fHistOS",100,massDplus-0.2,massDplus+0.2);
194 fHistOS->Sumw2();
195 fOutput->Add(fHistOS);
196
197 fHistOSbkg = new TH1F("fHistOSbkg","fHistOSbkg",100,massDplus-0.2,massDplus+0.2);
198 fHistOSbkg->Sumw2();
199 fOutput->Add(fHistOSbkg);
200
201 fHistLSnoweight = new TH1F("fHistLSnoweight","fHistLSnoweight",100,massDplus-0.2,massDplus+0.2);
202 fHistLSnoweight->Sumw2();
203 fOutput->Add(fHistLSnoweight);
204
205 fHistLSsinglecut = new TH1F("fHistLSsinglecut","fHistLSsinglecut",100,massDplus-0.2,massDplus+0.2);
206 fHistLSsinglecut->Sumw2();
207 fOutput->Add(fHistLSsinglecut);
208 }
209
1f4e9722 210 if(fFillNtuple){
211 OpenFile(2); // 2 is the slot number of the ntuple
4afc48a2 212 fNtupleDplus = new TNtuple("fNtupleDplus","D +","pdg:Px:Py:Pz:PtTrue:VxTrue:VyTrue:VzTrue:Ptpi:PtK:Ptpi2:PtRec:PointingAngle:DecLeng:VxRec:VyRec:VzRec:InvMass:sigvert:d0Pi:d0K:d0Pi2");
1f4e9722 213 }
214
d486095a 215 return;
216}
217
218//________________________________________________________________________
219void AliAnalysisTaskSEDplus::UserExec(Option_t */*option*/)
220{
221 // Execute analysis for current event:
222 // heavy flavor candidates association to MC truth
223
224
225
226 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
4afc48a2 227 fHistNEvents->Fill(0); // count event
228 // Post the data already here
229 PostData(1,fOutput);
230
d486095a 231
b557eb43 232 TClonesArray *array3Prong = 0;
4afc48a2 233 TClonesArray *arrayLikeSign =0;
b557eb43 234 if(!aod && AODEvent() && IsStandardAOD()) {
235 // In case there is an AOD handler writing a standard AOD, use the AOD
236 // event in memory rather than the input (ESD) event.
237 aod = dynamic_cast<AliAODEvent*> (AODEvent());
238 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
239 // have to taken from the AOD event hold by the AliAODExtension
240 AliAODHandler* aodHandler = (AliAODHandler*)
241 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
242 if(aodHandler->GetExtensions()) {
243 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
244 AliAODEvent *aodFromExt = ext->GetAOD();
245 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
4afc48a2 246 arrayLikeSign=(TClonesArray*)aodFromExt->GetList()->FindObject("LikeSign3Prong");
b557eb43 247 }
248 } else {
249 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
4afc48a2 250 arrayLikeSign=(TClonesArray*)aod->GetList()->FindObject("LikeSign3Prong");
b557eb43 251 }
8931c313 252
d486095a 253 if(!array3Prong) {
254 printf("AliAnalysisTaskSEDplus::UserExec: Charm3Prong branch not found!\n");
255 return;
256 }
4afc48a2 257 if(!arrayLikeSign) {
258 printf("AliAnalysisTaskSEDplus::UserExec: LikeSign3Prong branch not found!\n");
259 // do in any case the OS analysis.
260 }
261
262
263 TClonesArray *arrayMC=0;
264 AliAODMCHeader *mcHeader=0;
d486095a 265
266 // AOD primary vertex
1f4e9722 267 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
268 // vtx1->Print();
269
d486095a 270 // load MC particles
4afc48a2 271 if(fReadMC){
272
273 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
274 if(!arrayMC) {
275 printf("AliAnalysisTaskSEDplus::UserExec: MC particles branch not found!\n");
276 return;
277 }
1f4e9722 278
d486095a 279 // load MC header
4afc48a2 280 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
281 if(!mcHeader) {
282 printf("AliAnalysisTaskSEDplus::UserExec: MC header branch not found!\n");
283 return;
284 }
d486095a 285 }
4afc48a2 286
1f4e9722 287 Int_t n3Prong = array3Prong->GetEntriesFast();
288 if(fDebug>1) printf("Number of D+->Kpipi: %d\n",n3Prong);
4afc48a2 289
290
291 Int_t nOS=0;
1f4e9722 292 Int_t pdgDgDplustoKpipi[3]={321,211,211};
293 Double_t cutsDplus[12]={0.2,0.4,0.4,0.,0.,0.01,0.06,0.02,0.,0.85,0.,10000000000.};
294 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
295 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
296
297
298 Bool_t unsetvtx=kFALSE;
299 if(!d->GetOwnPrimaryVtx()){
300 d->SetOwnPrimaryVtx(vtx1);
301 unsetvtx=kTRUE;
302 }
4afc48a2 303
1f4e9722 304 if(d->SelectDplus(fVHF->GetDplusCuts())) {
305 Int_t iPtBin=0;
306 Double_t ptCand = d->Pt();
307 if(ptCand<2.){
308 iPtBin=0;
309 cutsDplus[7]=0.08;
310 cutsDplus[8]=0.5;
4afc48a2 311 cutsDplus[9]=0.979;
312 cutsDplus[10]=0.0000055;
1f4e9722 313 }
314 else if(ptCand>2. && ptCand<3){
315 iPtBin=1;
316 cutsDplus[7]=0.08;
317 cutsDplus[8]=0.5;
318 cutsDplus[9]=0.991;
4afc48a2 319 cutsDplus[10]=0.000005;
1f4e9722 320 }else if(ptCand>3. && ptCand<5){
321 iPtBin=2;
322 cutsDplus[7]=0.1;
323 cutsDplus[8]=0.5;
324 cutsDplus[9]=0.9955;
4afc48a2 325 cutsDplus[10]=0.0000035;
1f4e9722 326 }else{
327 iPtBin=3;
328 cutsDplus[7]=0.1;
329 cutsDplus[8]=0.5;
330 cutsDplus[9]=0.997;
4afc48a2 331 cutsDplus[10]=0.000001;
1f4e9722 332 }
333 Bool_t passTightCuts=d->SelectDplus(cutsDplus);
4afc48a2 334 Int_t labDp=-1;
335 Float_t deltaPx=0.;
336 Float_t deltaPy=0.;
337 Float_t deltaPz=0.;
338 Float_t truePt=0.;
339 Float_t xDecay=0.;
340 Float_t yDecay=0.;
341 Float_t zDecay=0.;
342 Float_t pdgCode=-2;
343 if(fReadMC){
344 labDp = d->MatchToMC(411,arrayMC,3,pdgDgDplustoKpipi);
345 if(labDp>=0){
346 AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
347 AliAODMCParticle *dg0 = (AliAODMCParticle*)arrayMC->At(partDp->GetDaughter(0));
348 deltaPx=partDp->Px()-d->Px();
349 deltaPy=partDp->Py()-d->Py();
350 deltaPz=partDp->Pz()-d->Pz();
351 truePt=partDp->Pt();
352 xDecay=dg0->Xv();
353 yDecay=dg0->Yv();
354 zDecay=dg0->Zv();
355 pdgCode=TMath::Abs(partDp->GetPdgCode());
356 }else{
357 pdgCode=-1;
358 }
359 }
1f4e9722 360 Double_t invMass=d->InvMassDplus();
361
362 TString hisNameA(Form("hMassPt%d",iPtBin));
363 TString hisNameS(Form("hSigPt%d",iPtBin));
364 TString hisNameB(Form("hBkgPt%d",iPtBin));
365 TString hisNameATC(Form("hMassPt%dTC",iPtBin));
366 TString hisNameSTC(Form("hSigPt%dTC",iPtBin));
367 TString hisNameBTC(Form("hBkgPt%dTC",iPtBin));
368
369 ((TH1F*)(fOutput->FindObject(hisNameA)))->Fill(invMass);
370 if(passTightCuts){
371 ((TH1F*)(fOutput->FindObject(hisNameATC)))->Fill(invMass);
372 }
d486095a 373
4afc48a2 374 Float_t tmp[22];
375 if(fFillNtuple){
376 tmp[0]=pdgCode;
377 tmp[1]=deltaPx;
378 tmp[2]=deltaPy;
379 tmp[3]=deltaPz;
380 tmp[4]=truePt;
381 tmp[5]=xDecay;
382 tmp[6]=yDecay;
383 tmp[7]=zDecay;
384 tmp[8]=d->PtProng(0);
385 tmp[9]=d->PtProng(1);
386 tmp[10]=d->PtProng(2);
387 tmp[11]=d->Pt();
388 tmp[12]=d->CosPointingAngle();
389 tmp[13]=d->DecayLength();
390 tmp[14]=d->Xv();
391 tmp[15]=d->Yv();
392 tmp[16]=d->Zv();
393 tmp[17]=d->InvMassDplus();
394 tmp[18]=d->GetSigmaVert();
395 tmp[19]=d->Getd0Prong(0);
396 tmp[20]=d->Getd0Prong(1);
397 tmp[21]=d->Getd0Prong(2);
398 fNtupleDplus->Fill(tmp);
399 PostData(2,fNtupleDplus);
400 }
d486095a 401
4afc48a2 402 if(fReadMC){
403 if(labDp>=0) {
404 ((TH1F*)(fOutput->FindObject(hisNameS)))->Fill(invMass);
405 if(passTightCuts){
406 ((TH1F*)(fOutput->FindObject(hisNameSTC)))->Fill(invMass);
407 }
408
409 }else{
410 ((TH1F*)(fOutput->FindObject(hisNameB)))->Fill(invMass);
411 if(passTightCuts){
412 ((TH1F*)(fOutput->FindObject(hisNameBTC)))->Fill(invMass);
413 }
414 fHistOSbkg->Fill(invMass);
fc8d975b 415 }
8c34bd86 416 }
4afc48a2 417
418 fHistOS->Fill(invMass);
419 nOS++;
d486095a 420 }
1f4e9722 421 if(unsetvtx) d->UnsetOwnPrimaryVtx();
422 }
4afc48a2 423
424 //start LS analysis
425 if(fDoLS && arrayLikeSign) LSAnalysis(array3Prong,arrayLikeSign,aod,vtx1,nOS);
426
427 PostData(1,fOutput);
d486095a 428 return;
429}
430
4afc48a2 431
432
433//_________________________________________________________________
434void AliAnalysisTaskSEDplus::LSAnalysis(TClonesArray *arrayOppositeSign,TClonesArray *arrayLikeSign,AliAODEvent *aod,AliAODVertex *vtx1, Int_t nDplusOS){
435
436/*
437 * Fill the Like Sign histograms
438 */
439
440
441
442
443 //count pos/neg tracks
444 Int_t nPosTrks=0,nNegTrks=0;
445 //counter for particles passing single particle cuts
446 Int_t nspcplus=0;
447 Int_t nspcminus=0;
448
449 for(Int_t it=0;it<aod->GetNumberOfTracks();it++) {
450 AliAODTrack *track = aod->GetTrack(it);
451 if(track->Charge()>0){
452 nPosTrks++;
453 if(track->Pt()>=0.4){
454 nspcplus++;
455 }
456 }else{
457 nNegTrks++;
458 if(track->Pt()>=0.4){
459 nspcminus++;
460 }
461 }
462 }
463
464 Int_t nOStriplets = arrayOppositeSign->GetEntriesFast();
465
466 Int_t nDplusLS=0;
467 Int_t nLikeSign = arrayLikeSign->GetEntriesFast();
468
469 for(Int_t iLikeSign = 0; iLikeSign < nLikeSign; iLikeSign++) {
470 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)arrayLikeSign->UncheckedAt(iLikeSign);
471 Bool_t unsetvtx=kFALSE;
472 if(!d->GetOwnPrimaryVtx()) {
473 d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
474 unsetvtx=kTRUE;
475 }
476 if(d->SelectDplus(fVHF->GetDplusCuts()))nDplusLS++;
477 if(unsetvtx) d->UnsetOwnPrimaryVtx();
478 }
479
480 Float_t wei2=0;
481 if(nLikeSign!=0)wei2 = (Float_t)nOStriplets/(Float_t)nLikeSign;
482 Float_t wei3=0;
483 if(nDplusLS!=0)wei3 = (Float_t)nDplusOS/(Float_t)nDplusLS;
484
485 // loop over sign candidates
486 for(Int_t iLikeSign = 0; iLikeSign < nLikeSign; iLikeSign++) {
487 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)arrayLikeSign->UncheckedAt(iLikeSign);
488 Bool_t unsetvtx=kFALSE;
489 if(!d->GetOwnPrimaryVtx()) {
490 d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
491 unsetvtx=kTRUE;
492 }
493
494 if(d->SelectDplus(fVHF->GetDplusCuts())){
495 Int_t sign= d->GetCharge();
496 Float_t wei1=1;
497 Float_t wei4=1;
498 if(sign>0&&nPosTrks>2&&nspcplus>2){//wei* should be automatically protected, since to get a triplet there must be at least 3 good tracks in the event
499 wei1=3.*(Float_t)nNegTrks/((Float_t)nPosTrks-2.);
500 wei4=3.*(Float_t)nspcminus/((Float_t)nspcplus-2.);
501 }
502 if(sign<0&&nNegTrks>2&&nspcminus>2){
503 wei1=3.*(Float_t)nPosTrks/((Float_t)nNegTrks-2.);
504 wei4=3.*(Float_t)nspcplus/((Float_t)nspcminus-2.);
505 }
506 Double_t invMass=d->InvMassDplus();
507 fHistLS->Fill(invMass,wei1);
508 fHistLSnoweight->Fill(invMass);
509 fHistLStrip->Fill(invMass,wei2);
510 fHistLStripcut->Fill(invMass,wei3);
511 fHistLSsinglecut->Fill(invMass,wei4);
512 }
513 if(unsetvtx) d->UnsetOwnPrimaryVtx();
514 }
515
516 //printf("------------ N. of positive tracks in Event ----- %d \n", nPosTrks);
517 //printf("------------ N. of negative tracks in Event ----- %d \n", nNegTrks);
518
519 // printf("LS analysis...done\n");
520
521}
522
523//_________________________________________________________________
524
525
d486095a 526void AliAnalysisTaskSEDplus::Terminate(Option_t */*option*/)
527{
528 // Terminate analysis
529 //
530 if(fDebug > 1) printf("AnalysisTaskSEDplus: Terminate() \n");
d486095a 531 fOutput = dynamic_cast<TList*> (GetOutputData(1));
532 if (!fOutput) {
533 printf("ERROR: fOutput not available\n");
534 return;
535 }
4afc48a2 536 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
1f4e9722 537 if(fFillNtuple){
538 fNtupleDplus = dynamic_cast<TNtuple*>(GetOutputData(2));
4afc48a2 539 }
540 fHistOS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistOS"));
541 fHistOSbkg = dynamic_cast<TH1F*>(fOutput->FindObject("fHistOSbkg"));
542 if(fDoLS){
543 fHistLS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistLS"));
544 fHistLStrip = dynamic_cast<TH1F*>(fOutput->FindObject("fHistLStrip"));
545 fHistLStripcut = dynamic_cast<TH1F*>(fOutput->FindObject("fHistLStripcut"));
546 fHistLSnoweight = dynamic_cast<TH1F*>(fOutput->FindObject("fHistLSnoweight"));
547 fHistLSsinglecut = dynamic_cast<TH1F*>(fOutput->FindObject("fHistLSsinglecut"));
1f4e9722 548 }
fc8d975b 549
4afc48a2 550
8c34bd86 551 return;
d486095a 552}