[u/mrichter/AliRoot.git] / FASTSIM / AliFastMuonTriggerEff.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */
// The trigger parametrization is computed for background levels 0., 0.5 and 1.
// In order to set a background level different from 0 it is necessary to 
// explicitly force it with:
// ForceBkgLevel(BkgLevel).
// For intermediate background levels, the trigger response is linearly 
// interpolated between these values.
// There is increased granularity in the pT region below 3 GeV. Although
// it does not seem to be necessary it is also possible to interpolate
// between pT bins using SetInt().


#include "AliFastMuonTriggerEff.h"
#include "TROOT.h"
#include "TFile.h"
#include "stdlib.h"
#include "TH3.h"
#include "TObjString.h"

#define PLIN printf("%s: %d: ",__FILE__,__LINE__)

ClassImp(AliFastMuonTriggerEff)

AliFastMuonTriggerEff::AliFastMuonTriggerEff():
    AliFastResponse("Efficiency", "Muon Trigger Efficiency")
{
//
// Default constructor
//
    SetCut(kLow);
    fZones=0;
    SetBkgLevel(0.);
    UnsetInt();
}

AliFastMuonTriggerEff::AliFastMuonTriggerEff(const AliFastMuonTriggerEff& eff)
    :AliFastResponse(eff)
{
// Copy constructor
    eff.Copy(*this);
}

void AliFastMuonTriggerEff::SetCut(Int_t cut) 
{  
  //
  // Set the pt cut
  if(cut==kLow){
    printf("Selecting Low Pt cut\n");
  }else if(cut==kHigh){
    printf("Selecting High Pt cut\n");
  }else if(cut==kAny){
    printf("Selecting Lowest Pt cut\n");
  }else{
    printf("Don't know cut %d! Selecting Low Pt cut\n",cut);
    cut=kLow;
  }
  fCut = cut;
  
}

Int_t AliFastMuonTriggerEff::SetBkgLevel(Float_t Bkg)
{
  //
  // Set the background level
  //
    if((Bkg!=0.)) {
       printf("%s: Warning: requested Bkg: %f\n",
       __FILE__,Bkg);
       fBkg=0.;
       printf("A consistent treatement of the trigger probability\n");
       printf("within the framework of the fast simulation requires\n");
       printf("requires background 0\n");
       printf("%s: fBkg: set to %f\n",
       __FILE__,fBkg);
    } else {
      fBkg=Bkg;
    }
    if(fZones!=0.) {
	Init();
    }
    return 0;
}

Int_t AliFastMuonTriggerEff::ForceBkgLevel(Float_t Bkg)
{
  //
  // Check and enforce consistency of the background level 
  // 
    if((Bkg!=0.)) {
       printf("%s: Warning: requested Bkg: %f\n",
       __FILE__,Bkg);
       printf("A consistent treatement of the trigger probability\n");
       printf("within the framework of the fast simulation\n");
       printf("requires background 0");
       printf("%s: Continue with fBkg: %f\n",
       __FILE__,Bkg);
    }
    fBkg=Bkg;
    if(fZones!=0.) {
	Init();
    }
    return 0;
}

Int_t AliFastMuonTriggerEff::LoadTables(Char_t *namet=""){
  //
  // Load the trigger tables
  //
    Char_t hNameA[100],hNameL[100],hNameH[100];
    sprintf(hNameA,"hEffAPt%s",namet);
    sprintf(hNameL,"hEffLPt%s",namet);
    sprintf(hNameH,"hEffHPt%s",namet);
    fhEffAPt = (TH3F*)gDirectory->Get(hNameA);
    fhEffLPt = (TH3F*)gDirectory->Get(hNameL);
    fhEffHPt = (TH3F*)gDirectory->Get(hNameH);
    if(!fhEffAPt){
      PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameA);
      return -1;
    }
    if(!fhEffLPt){
      PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameL);
      return -2;
    }
    if(!fhEffHPt){
      PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameH);
      return -3;
    }
    return 0;
}

void AliFastMuonTriggerEff::Init()
{
//
//  Initialization
//
    fZones=0;
    Char_t file[100]="$(ALICE_ROOT)/FASTSIM/data/MUONtriggerLUT_V2.4nvdn.root";
    printf("Initializing %s / %s\n", fName.Data(), fTitle.Data());
    printf("using data from file: %s\n",file);
    printf("AliFastMuonTriggerEff: Initialization with background level: %f\n",fBkg);
    TFile *f = new TFile(file);
    if(f->IsZombie()) {
        PLIN; printf("Cannot open file: %s\n",file);
        return;
    }
    f->ls();
    Int_t intb=0;
    Char_t namet[10];
    if(TMath::Abs(fBkg)<0.00001){
      sprintf(namet,"00");
    }else if(TMath::Abs(fBkg-0.5)<0.00001){
      sprintf(namet,"05");
    }else if(TMath::Abs(fBkg-1.0)<0.00001){
      sprintf(namet,"10");
    }else{
      PLIN; printf("A table for Bkg level: %f does not exists\n",fBkg);
      intb=1;
    }
    if(intb){ // Interpolation between background levels
      PLIN; printf("Interpolating Bkg level: %f\n",fBkg);
      TH3F* ha1,*hl1,*hh1,*ha2,*hl2,*hh2,*ha0,*hl0,*hh0;
      Char_t name1[10],name2[10]; Float_t b1,b2;
      if(fBkg>0&&fBkg<0.5){
        sprintf(name1,"00");
        sprintf(name2,"05");
	b1=0.;
	b2=0.5;
      }else if(fBkg>0.5){
        sprintf(name1,"05");
        sprintf(name2,"10");
	b1=0.5;
	b2=1.0;
	if(fBkg>1.0){
	  for(Int_t i=0; i<10;i++){
	    PLIN; printf("WARNING!!!! You are extrapolating above background 1.0\n");
	  }
	}
      }else{
        PLIN; printf("Bkg level: %f is not supported\n",fBkg);
        return;
      }
      if(LoadTables(name1)){
        PLIN; printf("Error in loading trigger tables\n");
	return;
      }
      PLIN; printf("We use tables for %f and %f to interpolate %f Bkg level\n",b1,b2,fBkg);
      ha0=(TH3F*)fhEffAPt->Clone("hEffAPtXX"); ha0->Reset();
      hl0=(TH3F*)fhEffLPt->Clone("hEffLPtXX"); hl0->Reset();
      hh0=(TH3F*)fhEffHPt->Clone("hEffHPtXX"); hh0->Reset();
      ha1=fhEffAPt;
      hl1=fhEffLPt;
      hh1=fhEffHPt;
      if(LoadTables(name2)){
        PLIN; printf("Error in loading trigger tables\n");
	return;
      }
      ha2=fhEffAPt;
      hl2=fhEffLPt;
      hh2=fhEffHPt;
      fhEffAPt=ha0;
      fhEffLPt=hl0;
      fhEffHPt=hh0;
      Int_t nnx=ha0->GetNbinsX()+1;
      Int_t nny=ha0->GetNbinsY()+1;
      Int_t nnz=ha0->GetNbinsZ()+1;
      for(Int_t ix=0; ix<=nnx; ix++){
	for(Int_t iy=0; iy<=nny; iy++){
	  for(Int_t iz=0; iz<=nnz; iz++){
	    Double_t y1,y2; Float_t cont;
	    y1=ha1->GetBinContent(ix,iy,iz); y2=ha2->GetBinContent(ix,iy,iz);
	    cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	    fhEffAPt->SetBinContent(ix,iy,iz,cont);
	    y1=hl1->GetBinContent(ix,iy,iz); y2=hl2->GetBinContent(ix,iy,iz);
	    cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	    fhEffLPt->SetBinContent(ix,iy,iz,cont);
	    y1=hh1->GetBinContent(ix,iy,iz); y2=hh2->GetBinContent(ix,iy,iz);
	    cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	    fhEffHPt->SetBinContent(ix,iy,iz,cont);
	  }
	}
      }
    }else{ // Use tables computed for selected backgound levels
      printf("Loading tables for background level: %f\n",fBkg);
      if(LoadTables(namet)){
        PLIN; printf("Error in loading trigger tables\n");
	return;
      }
    }
    fhEffAPt->SetDirectory(0);
    fhEffLPt->SetDirectory(0);
    fhEffHPt->SetDirectory(0);
    fhLX=fhEffLPt->GetXaxis();
    fhLY=fhEffLPt->GetYaxis();
    fhLZ=fhEffLPt->GetZaxis();
//
//
    if(f->Get("Description"))
    {
      fDescription=((TObjString*)f->Get("Description"))->GetString();
      printf("%s\n",fDescription.Data());
    }

    fThetaMin = fhEffLPt->GetXaxis()->GetXmin();
    fThetaMax = fhEffLPt->GetXaxis()->GetXmax();
    fnthetab=fhEffLPt->GetNbinsX();
    fDtheta   = (fThetaMax-fThetaMin)/fnthetab;

    fPhiMin   = fhEffLPt->GetYaxis()->GetXmin();
    fPhiMax   = fhEffLPt->GetYaxis()->GetXmax();
    fnphib=fhEffLPt->GetNbinsY();
    fDphi     = (fPhiMax-fPhiMin)/fnphib;

    fPtMin=fhEffLPt->GetZaxis()->GetXmin();
    fPtMax=fhEffLPt->GetZaxis()->GetXmax();
    fnptb=fhEffLPt->GetNbinsZ();
    fDpt      = (fPtMax-fPtMin)/fnptb;

    printf("%4d bins in theta [%f:%f]\n",fnthetab,fThetaMin,fThetaMax);
    printf("%4d bins in phi [%f:%f]\n",fnphib,fPhiMin,fPhiMax);
    printf("%4d bins in pt [%f:%f]\n",fnptb,fPtMin,fPtMax);

    fZones=fnthetab*fnphib;

    f->Close();
    if(fInt==0) {
      printf("Interpolation of trigger efficiencies is off!\n");
    } else {
      printf("Interpolation of trigger efficiencies is on!\n");
    }
}

void AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,Float_t theta,
              Float_t phi, Float_t& effLow, Float_t& effHigh, Float_t& effAny)
{
    //
    //  Trigger efficiency for pt, theta, phi (low, high and "any" cut)
    //
#ifdef MYTRIGDEBUG
    printf("Evaluate(ch=%2.0f, pt=%10.6f, theta=%7.2f, phi=%8.2f ...)\n",charge,pt,theta,phi);
#endif
    effLow=0.;
    effHigh=0.;
    effAny=0;
    if(fZones==0) {
        printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
	return;
    }
    if(pt<0) {
        printf("Warning: pt: %f < 0. GeV/c\n",pt);
	return;	
    }

    Int_t iPt   = fhLZ->FindBin((Double_t)pt);
    if(iPt>fnptb)iPt=fnptb;
    Int_t iPhi  = Int_t((phi-fPhiMin)/fDphi);
    if(phi<fPhiMin)iPhi=iPhi-1;
    Int_t iTheta = fhLX->FindBin((Double_t)theta);
#ifdef MYTRIGDEBUG
    printf("Evaluate(ch=%2.0f, pt=%10.6f, theta=%7.2f, phi=%8.2f ...)\n",charge,pt,theta,phi);
    printf(" 0:%1d iPt iTheta iPhi: %d %d %d\n",fInt,iPt,iTheta,iPhi);
#endif
    iPhi=iPhi-2*fnphib*(iPhi/(2*fnphib));
#ifdef MYTRIGDEBUG
    printf(" 1:%1d iPhi converted to: %d for angle equivalence\n",fInt,iPhi);
#endif
    if(iPhi<0)iPhi=-iPhi-1;
    if(iPhi>(fnphib-1))iPhi=2*fnphib-1-iPhi;
#ifdef MYTRIGDEBUG
    printf(" 2:%1d iPhi converted to: %d for the symmetry of the spectrometer\n",fInt,iPhi);
#endif
    if(charge==1.){
    } else if(charge==-1.) {
    iPhi=fnphib-1-iPhi;
#ifdef MYTRIGDEBUG
    printf(" 3:%1d iPhi converted to: %d for the charge symmetry\n",fInt,iPhi);
#endif
    } else {
        printf("Warning: not understand charge: %f\n",charge);
        return;
    }
    if(iTheta<=0||iTheta>fnthetab) {
        printf("Warning: theta: %f outside acceptance\n",theta);
        return;
    }
    if(iPt<0) {
        printf("Warning: what do you mean with pt: %f <0?\n",pt);
	return;
    }
    iPhi++;
#ifdef MYTRIGDEBUG
    printf(" 4:%1d Getting: iTheta, iPhi, iPt: %d %d %d\n",
              fInt,iTheta,iPhi,iPt);
#endif
    effLow =fhEffLPt->GetBinContent(iTheta,iPhi,iPt);
    effHigh=fhEffHPt->GetBinContent(iTheta,iPhi,iPt);
    effAny =fhEffAPt->GetBinContent(iTheta,iPhi,iPt);
#ifdef MYTRIGDEBUG
    printf(" 4:%1d Result: charge, iTheta, iPhi, iPt: %f %d %d %d effLow: %f, effHigh: %f, effAny: %f\n",
              fInt,charge,iTheta,iPhi,iPt,effLow,effHigh,effAny);
#endif
        
    if(fInt==1) {
      Float_t angl,angh,anga;
      Float_t effLowp,effHighp,effAnyp;
      Float_t ptc=(iPt+0.5)*fDpt;  // The center of current bin
      #ifdef MYTRIGDEBUG
        printf(" 5:1 The center of current bin iPt: %d is: %f\n",iPt,ptc);
      #endif
      if(iPt==fnptb) {
        #ifdef MYTRIGDEBUG
	printf(" 6:1 No more points above! No interpolation is needed!\n");
        #endif
	return;        
      }else if(ptc==pt){
        #ifdef MYTRIGDEBUG
	  printf(" 6:1 No interpolation is needed!\n");
        #endif
	return;
      }else if(ptc>pt){
	// Looking for previous point
        if(iPt>1) {
          effLowp =fhEffLPt->GetBinContent(iTheta,iPhi,iPt-1);
          effHighp=fhEffHPt->GetBinContent(iTheta,iPhi,iPt-1);
          effAnyp =fhEffAPt->GetBinContent(iTheta,iPhi,iPt-1);
          #ifdef MYTRIGDEBUG
	  printf(" 7:1 A simple look to previous point: %d: %f %f\n",iPt-1,effLowp,effHighp);
          #endif
	} else {
	  effLowp=0.;
	  effHighp=0.;
          effAnyp=0;
          #ifdef MYTRIGDEBUG
          printf(" 8:1 result is: %f %f %f\n",effLowp,effHighp,effAnyp);
	  #endif	  
	}
        angl=(effLow-effLowp)/fDpt;
        angh=(effHigh-effHighp)/fDpt;
        anga=(effAny-effAnyp)/fDpt;
      }else{
	// Looking for next point
	if(iPt<fnptb) {
          effLowp =fhEffLPt->GetBinContent(iTheta,iPhi,iPt+1);
          effHighp=fhEffHPt->GetBinContent(iTheta,iPhi,iPt+1);
          effAnyp =fhEffAPt->GetBinContent(iTheta,iPhi,iPt+1);
          #ifdef MYTRIGDEBUG
	  printf(" 7:1 A simple look to next point: %d: %f %f %f\n",iPt-1,effLowp,effHighp,effAnyp);
          #endif
	} else {
	  effLowp=effLow;
	  effHighp=effHigh;
	  effAnyp=effAny;
          #ifdef MYTRIGDEBUG
          printf(" 8:1 result is: pt: %f %f %f\n",effLowp,effHighp,effAnyp);
	  #endif	  
        }
        angl=(effLowp-effLow)/fDpt;
        angh=(effHighp-effHigh)/fDpt;
        anga=(effAnyp-effAny)/fDpt;
      }
      effLow=effLow+angl*(pt-ptc);
      effHigh=effHigh+angh*(pt-ptc);
      effAny=effAny+anga*(pt-ptc);
      #ifdef MYTRIGDEBUG
      printf(" 9:1 the interpolation coefficients are: %f %f %f\n",angl,angh,anga);
      #endif
  }
  #ifdef MYTRIGDEBUG
  printf("10:%1d effLow, effHigh=%f %f %f\n",fInt,effLow,effHigh,effAny);
  #endif
  return;
}


Float_t AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,
                   Float_t theta, Float_t phi)
{
    //
    // Trigger efficiency for pt, theta, phi depending of fCut
    // 
    if(fZones==0) {
        printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
	return 0.;
    }
    Float_t eff;
    Float_t effLow, effHigh, effAny;
    
    Evaluate(charge,pt,theta,phi,effLow,effHigh,effAny);
    if (fCut == kLow) 
	eff  = effLow;
    else if (fCut == kHigh)
	eff  = effHigh;
    else if (fCut == kAny)
        eff  = effAny;
    else
        eff  = 0;

    return eff;
}

AliFastMuonTriggerEff& AliFastMuonTriggerEff::operator=(const  AliFastMuonTriggerEff& rhs)
{
// Assignment operator
    rhs.Copy(*this);
    return *this;
}
Commit	Line	Data
6255180c	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, 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
a42548b0	16	/* $Id$ */
	17	// The trigger parametrization is computed for background levels 0., 0.5 and 1.
	18	// In order to set a background level different from 0 it is necessary to
	19	// explicitly force it with:
	20	// ForceBkgLevel(BkgLevel).
	21	// For intermediate background levels, the trigger response is linearly
	22	// interpolated between these values.
	23	// There is increased granularity in the pT region below 3 GeV. Although
	24	// it does not seem to be necessary it is also possible to interpolate
	25	// between pT bins using SetInt().
6255180c	26
490243f3	27
490243f3	28	#include "AliFastMuonTriggerEff.h"
2724a683	29	#include "TROOT.h"
	30	#include "TFile.h"
	31	#include "stdlib.h"
	32	#include "TH3.h"
	33	#include "TObjString.h"
	34
490243f3	35	#define PLIN printf("%s: %d: ",__FILE__,__LINE__)
6255180c	36
ecda8707	37	ClassImp(AliFastMuonTriggerEff)
6255180c	38
	39	AliFastMuonTriggerEff::AliFastMuonTriggerEff():
	40	AliFastResponse("Efficiency", "Muon Trigger Efficiency")
	41	{
	42	//
	43	// Default constructor
	44	//
	45	SetCut(kLow);
	46	fZones=0;
490243f3	47	SetBkgLevel(0.);
490243f3	48	UnsetInt();
6255180c	49	}
6255180c	50
a42548b0	51	AliFastMuonTriggerEff::AliFastMuonTriggerEff(const AliFastMuonTriggerEff& eff)
	52	:AliFastResponse(eff)
	53	{
	54	// Copy constructor
	55	eff.Copy(*this);
	56	}
	57
490243f3	58	void AliFastMuonTriggerEff::SetCut(Int_t cut)
490243f3	59	{
2724a683	60	//
2724a683	61	// Set the pt cut
490243f3	62	if(cut==kLow){
	63	printf("Selecting Low Pt cut\n");
	64	}else if(cut==kHigh){
	65	printf("Selecting High Pt cut\n");
	66	}else if(cut==kAny){
	67	printf("Selecting Lowest Pt cut\n");
	68	}else{
	69	printf("Don't know cut %d! Selecting Low Pt cut\n",cut);
	70	cut=kLow;
	71	}
	72	fCut = cut;
	73
	74	}
6255180c	75
490243f3	76	Int_t AliFastMuonTriggerEff::SetBkgLevel(Float_t Bkg)
6255180c	77	{
2724a683	78	//
	79	// Set the background level
	80	//
490243f3	81	if((Bkg!=0.)) {
	82	printf("%s: Warning: requested Bkg: %f\n",
	83	__FILE__,Bkg);
	84	fBkg=0.;
	85	printf("A consistent treatement of the trigger probability\n");
	86	printf("within the framework of the fast simulation requires\n");
	87	printf("requires background 0\n");
	88	printf("%s: fBkg: set to %f\n",
	89	__FILE__,fBkg);
	90	} else {
	91	fBkg=Bkg;
	92	}
	93	if(fZones!=0.) {
	94	Init();
	95	}
	96	return 0;
	97	}
6255180c	98
490243f3	99	Int_t AliFastMuonTriggerEff::ForceBkgLevel(Float_t Bkg)
490243f3	100	{
2724a683	101	//
	102	// Check and enforce consistency of the background level
	103	//
490243f3	104	if((Bkg!=0.)) {
	105	printf("%s: Warning: requested Bkg: %f\n",
	106	__FILE__,Bkg);
	107	printf("A consistent treatement of the trigger probability\n");
	108	printf("within the framework of the fast simulation\n");
	109	printf("requires background 0");
	110	printf("%s: Continue with fBkg: %f\n",
	111	__FILE__,Bkg);
	112	}
	113	fBkg=Bkg;
	114	if(fZones!=0.) {
	115	Init();
	116	}
	117	return 0;
	118	}
6255180c	119
490243f3	120	Int_t AliFastMuonTriggerEff::LoadTables(Char_t *namet=""){
2724a683	121	//
	122	// Load the trigger tables
	123	//
490243f3	124	Char_t hNameA[100],hNameL[100],hNameH[100];
	125	sprintf(hNameA,"hEffAPt%s",namet);
	126	sprintf(hNameL,"hEffLPt%s",namet);
	127	sprintf(hNameH,"hEffHPt%s",namet);
	128	fhEffAPt = (TH3F*)gDirectory->Get(hNameA);
	129	fhEffLPt = (TH3F*)gDirectory->Get(hNameL);
	130	fhEffHPt = (TH3F*)gDirectory->Get(hNameH);
	131	if(!fhEffAPt){
	132	PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameA);
	133	return -1;
	134	}
	135	if(!fhEffLPt){
	136	PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameL);
	137	return -2;
	138	}
	139	if(!fhEffHPt){
	140	PLIN; printf("%s: histogram %s not found\n",__FILE__,hNameH);
	141	return -3;
	142	}
	143	return 0;
6255180c	144	}
6255180c	145
490243f3	146	void AliFastMuonTriggerEff::Init()
6255180c	147	{
6255180c	148	//
490243f3	149	// Initialization
6255180c	150	//
490243f3	151	fZones=0;
	152	Char_t file[100]="$(ALICE_ROOT)/FASTSIM/data/MUONtriggerLUT_V2.4nvdn.root";
	153	printf("Initializing %s / %s\n", fName.Data(), fTitle.Data());
	154	printf("using data from file: %s\n",file);
	155	printf("AliFastMuonTriggerEff: Initialization with background level: %f\n",fBkg);
6255180c	156	TFile *f = new TFile(file);
6255180c	157	if(f->IsZombie()) {
490243f3	158	PLIN; printf("Cannot open file: %s\n",file);
6255180c	159	return;
	160	}
	161	f->ls();
490243f3	162	Int_t intb=0;
	163	Char_t namet[10];
	164	if(TMath::Abs(fBkg)<0.00001){
	165	sprintf(namet,"00");
	166	}else if(TMath::Abs(fBkg-0.5)<0.00001){
	167	sprintf(namet,"05");
	168	}else if(TMath::Abs(fBkg-1.0)<0.00001){
	169	sprintf(namet,"10");
	170	}else{
	171	PLIN; printf("A table for Bkg level: %f does not exists\n",fBkg);
	172	intb=1;
6255180c	173	}
490243f3	174	if(intb){ // Interpolation between background levels
	175	PLIN; printf("Interpolating Bkg level: %f\n",fBkg);
	176	TH3F* ha1,hl1,hh1,ha2,hl2,hh2,ha0,hl0,hh0;
	177	Char_t name1[10],name2[10]; Float_t b1,b2;
	178	if(fBkg>0&&fBkg<0.5){
	179	sprintf(name1,"00");
	180	sprintf(name2,"05");
	181	b1=0.;
	182	b2=0.5;
	183	}else if(fBkg>0.5){
	184	sprintf(name1,"05");
	185	sprintf(name2,"10");
	186	b1=0.5;
	187	b2=1.0;
	188	if(fBkg>1.0){
	189	for(Int_t i=0; i<10;i++){
	190	PLIN; printf("WARNING!!!! You are extrapolating above background 1.0\n");
	191	}
6255180c	192	}
490243f3	193	}else{
	194	PLIN; printf("Bkg level: %f is not supported\n",fBkg);
	195	return;
	196	}
	197	if(LoadTables(name1)){
	198	PLIN; printf("Error in loading trigger tables\n");
	199	return;
	200	}
	201	PLIN; printf("We use tables for %f and %f to interpolate %f Bkg level\n",b1,b2,fBkg);
	202	ha0=(TH3F*)fhEffAPt->Clone("hEffAPtXX"); ha0->Reset();
	203	hl0=(TH3F*)fhEffLPt->Clone("hEffLPtXX"); hl0->Reset();
	204	hh0=(TH3F*)fhEffHPt->Clone("hEffHPtXX"); hh0->Reset();
	205	ha1=fhEffAPt;
	206	hl1=fhEffLPt;
	207	hh1=fhEffHPt;
	208	if(LoadTables(name2)){
	209	PLIN; printf("Error in loading trigger tables\n");
	210	return;
	211	}
	212	ha2=fhEffAPt;
	213	hl2=fhEffLPt;
	214	hh2=fhEffHPt;
	215	fhEffAPt=ha0;
	216	fhEffLPt=hl0;
	217	fhEffHPt=hh0;
	218	Int_t nnx=ha0->GetNbinsX()+1;
	219	Int_t nny=ha0->GetNbinsY()+1;
	220	Int_t nnz=ha0->GetNbinsZ()+1;
	221	for(Int_t ix=0; ix<=nnx; ix++){
	222	for(Int_t iy=0; iy<=nny; iy++){
	223	for(Int_t iz=0; iz<=nnz; iz++){
	224	Double_t y1,y2; Float_t cont;
	225	y1=ha1->GetBinContent(ix,iy,iz); y2=ha2->GetBinContent(ix,iy,iz);
	226	cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	227	fhEffAPt->SetBinContent(ix,iy,iz,cont);
	228	y1=hl1->GetBinContent(ix,iy,iz); y2=hl2->GetBinContent(ix,iy,iz);
	229	cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	230	fhEffLPt->SetBinContent(ix,iy,iz,cont);
	231	y1=hh1->GetBinContent(ix,iy,iz); y2=hh2->GetBinContent(ix,iy,iz);
	232	cont=Float_t(y1+(y2-y1)/(b2-b1)*(fBkg-b1)); if(cont>1)cont=1; if(cont<0)cont=0;
	233	fhEffHPt->SetBinContent(ix,iy,iz,cont);
	234	}
6255180c	235	}
490243f3	236	}
	237	}else{ // Use tables computed for selected backgound levels
	238	printf("Loading tables for background level: %f\n",fBkg);
	239	if(LoadTables(namet)){
	240	PLIN; printf("Error in loading trigger tables\n");
	241	return;
	242	}
6255180c	243	}
490243f3	244	fhEffAPt->SetDirectory(0);
	245	fhEffLPt->SetDirectory(0);
	246	fhEffHPt->SetDirectory(0);
	247	fhLX=fhEffLPt->GetXaxis();
	248	fhLY=fhEffLPt->GetYaxis();
	249	fhLZ=fhEffLPt->GetZaxis();
	250	//
	251	//
	252	if(f->Get("Description"))
	253	{
	254	fDescription=((TObjString*)f->Get("Description"))->GetString();
	255	printf("%s\n",fDescription.Data());
6255180c	256	}
490243f3	257
	258	fThetaMin = fhEffLPt->GetXaxis()->GetXmin();
	259	fThetaMax = fhEffLPt->GetXaxis()->GetXmax();
	260	fnthetab=fhEffLPt->GetNbinsX();
	261	fDtheta = (fThetaMax-fThetaMin)/fnthetab;
	262
	263	fPhiMin = fhEffLPt->GetYaxis()->GetXmin();
	264	fPhiMax = fhEffLPt->GetYaxis()->GetXmax();
	265	fnphib=fhEffLPt->GetNbinsY();
	266	fDphi = (fPhiMax-fPhiMin)/fnphib;
	267
	268	fPtMin=fhEffLPt->GetZaxis()->GetXmin();
	269	fPtMax=fhEffLPt->GetZaxis()->GetXmax();
	270	fnptb=fhEffLPt->GetNbinsZ();
	271	fDpt = (fPtMax-fPtMin)/fnptb;
	272
	273	printf("%4d bins in theta [%f:%f]\n",fnthetab,fThetaMin,fThetaMax);
	274	printf("%4d bins in phi [%f:%f]\n",fnphib,fPhiMin,fPhiMax);
	275	printf("%4d bins in pt [%f:%f]\n",fnptb,fPtMin,fPtMax);
	276
	277	fZones=fnthetab*fnphib;
	278
6255180c	279	f->Close();
490243f3	280	if(fInt==0) {
	281	printf("Interpolation of trigger efficiencies is off!\n");
	282	} else {
	283	printf("Interpolation of trigger efficiencies is on!\n");
	284	}
6255180c	285	}
6255180c	286
490243f3	287	void AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,Float_t theta,
490243f3	288	Float_t phi, Float_t& effLow, Float_t& effHigh, Float_t& effAny)
6255180c	289	{
490243f3	290	//
	291	// Trigger efficiency for pt, theta, phi (low, high and "any" cut)
	292	//
	293	#ifdef MYTRIGDEBUG
	294	printf("Evaluate(ch=%2.0f, pt=%10.6f, theta=%7.2f, phi=%8.2f ...)\n",charge,pt,theta,phi);
	295	#endif
6255180c	296	effLow=0.;
6255180c	297	effHigh=0.;
490243f3	298	effAny=0;
6255180c	299	if(fZones==0) {
	300	printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
	301	return;
	302	}
	303	if(pt<0) {
	304	printf("Warning: pt: %f < 0. GeV/c\n",pt);
	305	return;
	306	}
490243f3	307
	308	Int_t iPt = fhLZ->FindBin((Double_t)pt);
	309	if(iPt>fnptb)iPt=fnptb;
6255180c	310	Int_t iPhi = Int_t((phi-fPhiMin)/fDphi);
6255180c	311	if(phi<fPhiMin)iPhi=iPhi-1;
490243f3	312	Int_t iTheta = fhLX->FindBin((Double_t)theta);
6255180c	313	#ifdef MYTRIGDEBUG
490243f3	314	printf("Evaluate(ch=%2.0f, pt=%10.6f, theta=%7.2f, phi=%8.2f ...)\n",charge,pt,theta,phi);
490243f3	315	printf(" 0:%1d iPt iTheta iPhi: %d %d %d\n",fInt,iPt,iTheta,iPhi);
6255180c	316	#endif
490243f3	317	iPhi=iPhi-2fnphib(iPhi/(2*fnphib));
6255180c	318	#ifdef MYTRIGDEBUG
490243f3	319	printf(" 1:%1d iPhi converted to: %d for angle equivalence\n",fInt,iPhi);
6255180c	320	#endif
6255180c	321	if(iPhi<0)iPhi=-iPhi-1;
490243f3	322	if(iPhi>(fnphib-1))iPhi=2*fnphib-1-iPhi;
6255180c	323	#ifdef MYTRIGDEBUG
490243f3	324	printf(" 2:%1d iPhi converted to: %d for the symmetry of the spectrometer\n",fInt,iPhi);
6255180c	325	#endif
	326	if(charge==1.){
	327	} else if(charge==-1.) {
490243f3	328	iPhi=fnphib-1-iPhi;
6255180c	329	#ifdef MYTRIGDEBUG
490243f3	330	printf(" 3:%1d iPhi converted to: %d for the charge symmetry\n",fInt,iPhi);
6255180c	331	#endif
	332	} else {
	333	printf("Warning: not understand charge: %f\n",charge);
	334	return;
	335	}
490243f3	336	if(iTheta<=0\|\|iTheta>fnthetab) {
6255180c	337	printf("Warning: theta: %f outside acceptance\n",theta);
	338	return;
	339	}
	340	if(iPt<0) {
	341	printf("Warning: what do you mean with pt: %f <0?\n",pt);
	342	return;
	343	}
490243f3	344	iPhi++;
6255180c	345	#ifdef MYTRIGDEBUG
490243f3	346	printf(" 4:%1d Getting: iTheta, iPhi, iPt: %d %d %d\n",
490243f3	347	fInt,iTheta,iPhi,iPt);
6255180c	348	#endif
490243f3	349	effLow =fhEffLPt->GetBinContent(iTheta,iPhi,iPt);
	350	effHigh=fhEffHPt->GetBinContent(iTheta,iPhi,iPt);
	351	effAny =fhEffAPt->GetBinContent(iTheta,iPhi,iPt);
6255180c	352	#ifdef MYTRIGDEBUG
490243f3	353	printf(" 4:%1d Result: charge, iTheta, iPhi, iPt: %f %d %d %d effLow: %f, effHigh: %f, effAny: %f\n",
490243f3	354	fInt,charge,iTheta,iPhi,iPt,effLow,effHigh,effAny);
6255180c	355	#endif
490243f3	356
	357	if(fInt==1) {
	358	Float_t angl,angh,anga;
	359	Float_t effLowp,effHighp,effAnyp;
	360	Float_t ptc=(iPt+0.5)*fDpt; // The center of current bin
	361	#ifdef MYTRIGDEBUG
	362	printf(" 5:1 The center of current bin iPt: %d is: %f\n",iPt,ptc);
	363	#endif
	364	if(iPt==fnptb) {
	365	#ifdef MYTRIGDEBUG
	366	printf(" 6:1 No more points above! No interpolation is needed!\n");
	367	#endif
	368	return;
	369	}else if(ptc==pt){
	370	#ifdef MYTRIGDEBUG
	371	printf(" 6:1 No interpolation is needed!\n");
	372	#endif
	373	return;
	374	}else if(ptc>pt){
	375	// Looking for previous point
	376	if(iPt>1) {
	377	effLowp =fhEffLPt->GetBinContent(iTheta,iPhi,iPt-1);
	378	effHighp=fhEffHPt->GetBinContent(iTheta,iPhi,iPt-1);
	379	effAnyp =fhEffAPt->GetBinContent(iTheta,iPhi,iPt-1);
	380	#ifdef MYTRIGDEBUG
	381	printf(" 7:1 A simple look to previous point: %d: %f %f\n",iPt-1,effLowp,effHighp);
	382	#endif
	383	} else {
	384	effLowp=0.;
	385	effHighp=0.;
	386	effAnyp=0;
	387	#ifdef MYTRIGDEBUG
	388	printf(" 8:1 result is: %f %f %f\n",effLowp,effHighp,effAnyp);
	389	#endif
	390	}
	391	angl=(effLow-effLowp)/fDpt;
	392	angh=(effHigh-effHighp)/fDpt;
	393	anga=(effAny-effAnyp)/fDpt;
	394	}else{
	395	// Looking for next point
	396	if(iPt<fnptb) {
	397	effLowp =fhEffLPt->GetBinContent(iTheta,iPhi,iPt+1);
	398	effHighp=fhEffHPt->GetBinContent(iTheta,iPhi,iPt+1);
	399	effAnyp =fhEffAPt->GetBinContent(iTheta,iPhi,iPt+1);
	400	#ifdef MYTRIGDEBUG
	401	printf(" 7:1 A simple look to next point: %d: %f %f %f\n",iPt-1,effLowp,effHighp,effAnyp);
	402	#endif
	403	} else {
	404	effLowp=effLow;
	405	effHighp=effHigh;
	406	effAnyp=effAny;
	407	#ifdef MYTRIGDEBUG
	408	printf(" 8:1 result is: pt: %f %f %f\n",effLowp,effHighp,effAnyp);
	409	#endif
	410	}
	411	angl=(effLowp-effLow)/fDpt;
	412	angh=(effHighp-effHigh)/fDpt;
	413	anga=(effAnyp-effAny)/fDpt;
	414	}
	415	effLow=effLow+angl*(pt-ptc);
	416	effHigh=effHigh+angh*(pt-ptc);
	417	effAny=effAny+anga*(pt-ptc);
	418	#ifdef MYTRIGDEBUG
	419	printf(" 9:1 the interpolation coefficients are: %f %f %f\n",angl,angh,anga);
420	#endif
421	}
422	#ifdef MYTRIGDEBUG
423	printf("10:%1d effLow, effHigh=%f %f %f\n",fInt,effLow,effHigh,effAny);
424	#endif
425	return;
426	}
427
bc77d57e	428
6255180c	429
	430	Float_t AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,
	431	Float_t theta, Float_t phi)
	432	{
	433	//
	434	// Trigger efficiency for pt, theta, phi depending of fCut
	435	//
	436	if(fZones==0) {
	437	printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
	438	return 0.;
	439	}
	440	Float_t eff;
490243f3	441	Float_t effLow, effHigh, effAny;
6255180c	442
490243f3	443	Evaluate(charge,pt,theta,phi,effLow,effHigh,effAny);
6255180c	444	if (fCut == kLow)
6255180c	445	eff = effLow;
490243f3	446	else if (fCut == kHigh)
6255180c	447	eff = effHigh;
490243f3	448	else if (fCut == kAny)
	449	eff = effAny;
	450	else
	451	eff = 0;
6255180c	452
	453	return eff;
	454	}
a42548b0	455
	456	AliFastMuonTriggerEff& AliFastMuonTriggerEff::operator=(const AliFastMuonTriggerEff& rhs)
	457	{
	458	// Assignment operator
	459	rhs.Copy(*this);
	460	return *this;
	461	}
	462