SFHMarana: TTriggerReweight.C Source File

00001 
00002 // Class TTriggerReweight
00003 //
00004 // Author     : E. Sauvan
00005 // Created    : 24/11/2004
00006 // Last update: $Date: 2005/06/08 11:39:57 $
00007 //          by: $Author: blist $
00008 // Comment: Class to handle reweigthing from trigger efficiencies
00010 
00011 #include "TTriggerReweight.h"
00012 
00013 TTriggerReweight::TTriggerReweight()
00014 {
00015 
00016 fCCTrigHERA2Init=kTRUE;
00017 
00018 };
00019 
00020 
00021 
00022 Float_t TTriggerReweight::DiJetsTriggerEff(Float_t PtJet1,Float_t ThetaJet1,Int_t RunPeriod)
00023 {
00024 //--- Trigger effi for jet trigger (2-jet events) as a function of
00025 //      Pt/Theta of leading jet, run period
00026 //      From generic analysis (By Matti Peez)
00027   
00028   // Number of bins in pt and theta
00029   const Int_t NBinsPt = 5;
00030   const Int_t NBinsTh = 3;
00031   const Float_t ArrBinPt[NBinsPt] = {20.0,25.0,30.0,35.0,40.0};
00032   const Float_t ArrBinTh[NBinsTh] = {10.0,20.0,40.0};
00033   const Float_t MYR2D=180.0/3.14159265359;
00034   
00035   // Table of trigger efficiencies
00036   // Pt bins: 20,25,30,35,40,inf
00037   // Th bins: 10,20,40,150
00038   Float_t TriggEffBinThBinPt[3][5] = 
00039     {{1.0, 1.0, 1.0, 1.0, 1.0},
00040      {0.987, 1.0, 1.0, 1.0, 1.0},
00041      {0.986667, 1.0, 1.0 , 1.0 , 1.0}
00042     };
00043 
00044 //---- fill constantes for each run period 
00045 //     Y9497 = 3, Y9899 = 4 , Y9900 = 5, Y0203 = 6, Y0304 = 7
00046   if(RunPeriod==5) {//---- 9900
00047         TriggEffBinThBinPt[0][0]=1.0;
00048         TriggEffBinThBinPt[0][1]=1.0;
00049         TriggEffBinThBinPt[0][2]=1.0;
00050         TriggEffBinThBinPt[0][3]=1.0;
00051         TriggEffBinThBinPt[0][4]=1.0;
00052         TriggEffBinThBinPt[1][0]=0.987;
00053         TriggEffBinThBinPt[1][1]=1.0;
00054         TriggEffBinThBinPt[1][2]=1.0;
00055         TriggEffBinThBinPt[1][3]=1.0;
00056         TriggEffBinThBinPt[1][4]=1.0;
00057         TriggEffBinThBinPt[2][0]=0.98667;
00058         TriggEffBinThBinPt[3][1]=1.0;
00059         TriggEffBinThBinPt[4][2]=1.0;
00060         TriggEffBinThBinPt[5][3]=1.0;
00061         TriggEffBinThBinPt[6][4]=1.0;
00062   }else if(RunPeriod==3) {//---- 9407
00063         TriggEffBinThBinPt[0][0]=0.926829;
00064         TriggEffBinThBinPt[0][1]=0.926829;
00065         TriggEffBinThBinPt[0][2]=1.0;
00066         TriggEffBinThBinPt[0][3]=1.0;
00067         TriggEffBinThBinPt[0][4]=1.0;
00068         TriggEffBinThBinPt[1][0]=0.983051;
00069         TriggEffBinThBinPt[1][1]=0.983051;
00070         TriggEffBinThBinPt[1][2]=1.0;
00071         TriggEffBinThBinPt[1][3]=1.0;
00072         TriggEffBinThBinPt[1][4]=1.0;
00073         TriggEffBinThBinPt[2][0]=0.865672;
00074         TriggEffBinThBinPt[3][1]=0.958904;
00075         TriggEffBinThBinPt[4][2]=1.0;
00076         TriggEffBinThBinPt[5][3]=1.0;
00077         TriggEffBinThBinPt[6][4]=1.0;
00078   }else if(RunPeriod==4) {//---- 9899
00079         TriggEffBinThBinPt[0][0]=0.975;
00080         TriggEffBinThBinPt[0][1]=0.98;
00081         TriggEffBinThBinPt[0][2]=1.0;
00082         TriggEffBinThBinPt[0][3]=1.0;
00083         TriggEffBinThBinPt[0][4]=1.0;
00084         TriggEffBinThBinPt[1][0]=0.975;
00085         TriggEffBinThBinPt[1][1]=0.98;
00086         TriggEffBinThBinPt[1][2]=1.0;
00087         TriggEffBinThBinPt[1][3]=1.0;
00088         TriggEffBinThBinPt[1][4]=1.0;
00089         TriggEffBinThBinPt[2][0]=0.975;
00090         TriggEffBinThBinPt[3][1]=0.98;
00091         TriggEffBinThBinPt[4][2]=1.0;
00092         TriggEffBinThBinPt[5][3]=1.0;
00093         TriggEffBinThBinPt[6][4]=1.0;
00094   }else{//---- default set everything to 1
00095         TriggEffBinThBinPt[0][0]=1.0;
00096         TriggEffBinThBinPt[0][1]=1.0;
00097         TriggEffBinThBinPt[0][2]=1.0;
00098         TriggEffBinThBinPt[0][3]=1.0;
00099         TriggEffBinThBinPt[0][4]=1.0;
00100         TriggEffBinThBinPt[1][0]=0.987;
00101         TriggEffBinThBinPt[1][1]=1.0;
00102         TriggEffBinThBinPt[1][2]=1.0;
00103         TriggEffBinThBinPt[1][3]=1.0;
00104         TriggEffBinThBinPt[1][4]=1.0;
00105         TriggEffBinThBinPt[2][0]=0.98667;
00106         TriggEffBinThBinPt[3][1]=1.0;
00107         TriggEffBinThBinPt[4][2]=1.0;
00108         TriggEffBinThBinPt[5][3]=1.0;
00109         TriggEffBinThBinPt[6][4]=1.0;
00110   }
00111   
00112   
00113   // Final trigger efficiency depending of pt/th of highest pt jet
00114   Float_t TriggEff = 1.0;
00115 
00116 
00117   // Get trigger efficiency out of bank
00118   //  printf("TRIGGER: PT = %f TH = %f\n",PtJet1,ThetaJet1*MYR2D);
00119   Int_t PtBin = -1;
00120   Int_t ThBin = -1;
00121   // Find corresponding bin in pt to PtJet1
00122   for (Int_t iPt=0;iPt<NBinsPt;iPt++){
00123     if(iPt<NBinsPt-1){
00124       if(PtJet1>ArrBinPt[iPt] && PtJet1<ArrBinPt[iPt+1]){
00125         PtBin = iPt;
00126         break;
00127       }
00128     }else{ // Last bin => special treatment
00129       if(PtJet1>ArrBinPt[iPt]){
00130         PtBin = iPt;
00131         break;
00132       }
00133     }                  
00134   }
00135 
00136   // Find corresponding bin in theta to ThetaJet1
00137   for (Int_t iTh=0;iTh<NBinsTh;iTh++){
00138     if(iTh<NBinsTh-1){
00139       if(ThetaJet1*MYR2D>ArrBinTh[iTh] && ThetaJet1*MYR2D<ArrBinTh[iTh+1]){
00140         ThBin = iTh;
00141         break;
00142       }
00143     }else{// Last bin => special treatment
00144       if(ThetaJet1*MYR2D>ArrBinTh[iTh]){
00145         ThBin = iTh;
00146         break;
00147       }
00148     }                  
00149   }
00150 
00151   // Pt bins: 20,25,30,35,40,inf
00152   // Th bins: 10,20,40,150
00153   if((PtBin >= 0 && PtBin < NBinsPt) && (ThBin >= 0 && ThBin < NBinsTh)){
00154     TriggEff = TriggEffBinThBinPt[ThBin][PtBin];
00155   }
00156   // Else return 1.0
00157   else{
00158     return(1.0);
00159   }
00160 
00161   // Return trigger efficiency
00162   //  printf("Element Th = %d Pt = %d Eff = %f\n",ThBin,PtBin,TriggEff);
00163   return(TriggEff);
00164 
00165 
00166 }
00167 
00168 Float_t TTriggerReweight::CCTriggerEff(Float_t Pth,Float_t Thh,Int_t RunPeriod)
00169 {
00170 //---- Calc correction factors for CC trigger inefficiency
00171  const Float_t MYR2D=180.0/3.14159265359;
00172 
00173  if(RunPeriod==3){ //---9497: no trigger eff. available ..
00174         return(1.0);
00175  
00176  }else if(RunPeriod==4){ //---9899
00177 // Generated with f2c on trig_corr.F from W analysis (M. Schneider)
00178 // -> calls trig_corr__ converted by f2c
00179         Float_t TrigEff = 1.0;
00180 
00181         // Get trigger efficiency
00182         // x = PtH; y = ThH; imeth = 4;iyear = 98;
00183         Float_t PtHadrons = Pth;
00184         Float_t ThHadrons = Thh*MYR2D;
00185         Int_t IMET = 4;
00186         Int_t YEAR = 98;
00187         TrigEff = this->trigeff_9899__(&PtHadrons,&ThHadrons,&IMET,&YEAR);
00188 
00189         // Return trigger efficiency
00190         return(TrigEff);
00191  }else if(RunPeriod==5){ //---9900 
00192 // Generated with f2c on trig_corr.F from W analysis (M. Schneider)
00193 // -> calls trig_corr__ converted by f2c
00194         Float_t TrigEff = 1.0;
00195 
00196         // Init arrays (gammah bins)
00197         tr99_1.gbins[1] = 0;tr99_1.gbins[1] = 10.00;tr99_1.gbins[2] = 20.0;
00198         tr99_1.gbins[3] = 30.0;tr99_1.gbins[4] = 50.0;tr99_1.gbins[5] = 70.0;
00199         tr99_1.gbins[6] = 90.0;tr99_1.gbins[7] = 110.0;tr99_1.gbins[8] = 180.0;
00200 
00201         // constants
00202         tr99_1.out3[0] = 7.364261;tr99_1.out3[1] = 7.473708;tr99_1.out3[2] = 0.002017514;
00203         tr99_1.out3[3] = 0.9444522;tr99_1.out3[4] = 2.37768;tr99_1.out3[5] = 0.1343088;
00204         tr99_1.out3[6] = 1.003783;tr99_1.out3[7] = 3.14458;tr99_1.out3[8] = 0.1911239;
00205         tr99_1.out3[9] = 0.9975498;tr99_1.out3[10] = 3.832504;tr99_1.out3[11] = 0.2087284;
00206         tr99_1.out3[12] = 0.9980307;tr99_1.out3[13] = 2.966897;tr99_1.out3[14] = 0.1585862;
00207         tr99_1.out3[15] = 0.9715002;tr99_1.out3[16] = 2.252059;tr99_1.out3[17] = 0.1269788;
00208         tr99_1.out3[18] = 1.092006;tr99_1.out3[19] = 1.515356;tr99_1.out3[20] = 0.0648907;
00209         tr99_1.out3[21] = 1.22744;tr99_1.out3[22] = 1.321669;tr99_1.out3[23] = 0.04720018;
00210 
00211         // Get trigger efficiency
00212         // x = PtH; y = ThH; imeth = 4;iyear = 98;
00213         Float_t PtHadrons = Pth;
00214         Float_t ThHadrons = Thh*MYR2D;
00215         Int_t YEAR = 99;
00216         TrigEff = this->trig_corr__(&PtHadrons,&ThHadrons,&YEAR);
00217 
00218         // Return trigger efficiency
00219         return(TrigEff);
00220  }else return(1.0); //--- default: nothing ...
00221  
00222  
00223 }
00224 
00225 
00226 
00227 Float_t TTriggerReweight::ApplyCCTRIGGEFFreweightHERA2(Float_t xhadrons, Float_t Q2hadrons)
00228 {
00229 //      CC trigger efficiencies for 2004, By Benji...
00230 // 
00231 
00232   Float_t TriggCCWeight=1.0;
00233   
00234   if(fCCTrigHERA2Init==kTRUE){
00235         SetTriggEffTable();
00236         fCCTrigHERA2Init=kFALSE;
00237   }
00238                 
00239     //  Float_t efficiency=GetPsCCTriggerEfficiency(event.hfs);
00240     if(xhadrons<=0 || Q2hadrons<=0){
00241       TriggCCWeight=1;
00242       return(TriggCCWeight);
00243     }
00244 
00245 
00246  Double_t logxh=TMath::Log10(xhadrons);
00247  Double_t logQ2h=TMath::Log10(Q2hadrons);
00248 
00249     Bool_t foundx=kFALSE;
00250     Bool_t foundq2=kFALSE;
00251 
00252     Int_t ibinx=7;
00253     Int_t ibinq2=7;
00254 
00255  for(Int_t i=0;i<8;i++){
00256    if( xl[i] < logxh &&  logxh <= xh[i]  ) {
00257      ibinx=i;
00258      foundx=kTRUE;
00259    }
00260  if( q2l[i] < logQ2h &&  logQ2h <= q2h[i]  ) {
00261      ibinq2=i;
00262      foundq2=kTRUE;
00263    }
00264  }
00265 
00266  // efficiency=triggconst[ibinx][inbinq2];
00267  
00268  if (foundx==kFALSE) {
00269    if (logxh <= xl[0] ) ibinx=0;
00270    if (logxh > xh[7] ) ibinx=7;
00271  }
00272  
00273  if (foundq2==kFALSE) {
00274    if (logQ2h <= q2l[0] ) ibinq2=0;
00275    if (logQ2h > q2h[7] ) ibinq2=7;
00276  }
00277 
00278 
00279 
00280 //  for(Int_t i=0;i<8;i++){
00281 
00282 //    cout <<triggconst[i][0] <<" "<< triggconst[i][1] <<" " << triggconst[i][2] <<" " << triggconst[i][3]<< " " << triggconst[i][4] << " " << triggconst[i][5]<< " " << triggconst[i][6] << " " << triggconst[i][7]  << endl;
00283 
00284 // }
00285 
00286  Double_t efficiency=triggconst[ibinx][ibinq2];
00287 
00288  //   cout << "ibinx " << ibinx << " ibinQ2 " << ibinq2 << endl;
00289  //  cout << "efficiency " << efficiency <<    " log10(xh) " << logxh <<" log10(Q2h) " << logQ2h << endl;
00290 
00291   TriggCCWeight=efficiency;
00292 
00293 // //--- For systematics calculation: error on the efficiency 
00294 //   if(event.Syst_Type==13) {
00295 //         Float_t EffError=0.15*(1.-efficiency)+0.02; 
00296 //      event.TriggCCWeight=efficiency*(1+event.Syst_Sign*EffError);
00297 //   }
00298 
00299   return TriggCCWeight;
00300 }
00301 
00302 
00303 void TTriggerReweight::SetTriggEffTable()
00304 {
00305 
00306   xl[0]= -2.33; 
00307   xl[1]= -2.0; 
00308   xl[2]=  -1.67; 
00309   xl[3]=  -1.33; 
00310   xl[4]=  -1.0 ; 
00311   xl[5]=  -0.75 ; 
00312   xl[6]=  -0.5; 
00313   xl[7]= -0.25 ; 
00314  
00315   xh[0]= -2.0;
00316   xh[1]= -1.67;
00317   xh[2]= -1.33;
00318   xh[3]= -1.0;
00319   xh[4]= -0.75;
00320   xh[5]=  -0.5;
00321   xh[6]=  -0.25;
00322   xh[7]=  -0.0;
00323 
00324   q2l[0]= 2.35;
00325   q2l[1]=    2.6;
00326   q2l[2]=  2.85;
00327   q2l[3]=  3.1;
00328   q2l[4]=  3.35;
00329   q2l[5]= 3.6;
00330   q2l[6]= 3.85;
00331   q2l[7]= 4.1 ;
00332     
00333   q2h[0]=    2.6;
00334   q2h[1]=  2.85;
00335   q2h[2]=  3.1;
00336   q2h[3]=  3.35;
00337   q2h[4]= 3.6;
00338   q2h[5]= 3.85;
00339   q2h[6]= 4.1 ;
00340   q2h[7]= 4.4 ;
00341 
00342 
00343 
00344 triggconst[0][0]=0.356411;
00345 triggconst[0][1]=0.555255;
00346 triggconst[0][2]=1;
00347 triggconst[0][3]=1;
00348 triggconst[0][4]=1;
00349 triggconst[0][5]=1;
00350 triggconst[0][6]=1;
00351 triggconst[0][7]=1; 
00352 
00353 triggconst[1][0]=0.59701;
00354 triggconst[1][1]=0.686556;
00355 triggconst[1][2]=0.658748;
00356 triggconst[1][3]=0.778839;
00357 triggconst[1][4]=1;
00358 triggconst[1][5]=1;
00359 triggconst[1][6]=1;
00360 triggconst[1][7]=1; 
00361 
00362 triggconst[2][0]=0.776412;
00363 triggconst[2][1]=0.883318;
00364 triggconst[2][2]=0.930311;
00365 triggconst[2][3]=0.896854;
00366 triggconst[2][4]=0.917348;
00367 triggconst[2][5]=1;
00368 triggconst[2][6]=1;
00369 triggconst[2][7]=1; 
00370 
00371 triggconst[3][0]=0.887515;
00372 triggconst[3][1]=0.947987;
00373 triggconst[3][2]=0.980887;
00374 triggconst[3][3]=0.994757;
00375 triggconst[3][4]=0.987839;
00376 triggconst[3][5]=0.9588;
00377 triggconst[3][6]=1;
00378 triggconst[3][7]=1; 
00379 
00380 triggconst[4][0]=0.750268;
00381 triggconst[4][1]=0.816743;
00382 triggconst[4][2]=0.932612;
00383 triggconst[4][3]=1;
00384 triggconst[4][4]=1;
00385 triggconst[4][5]=1;
00386 triggconst[4][6]=1;
00387 triggconst[4][7]=1; 
00388 
00389 triggconst[5][0]=1;
00390 triggconst[5][1]=0.166138;
00391 triggconst[5][2]=0.556489;
00392 triggconst[5][3]=0.949651;
00393 triggconst[5][4]=1;
00394 triggconst[5][5]=1;
00395 triggconst[5][6]=1;
00396 triggconst[5][7]=1; 
00397 
00398 triggconst[6][0]=1;
00399 triggconst[6][1]=1;
00400 triggconst[6][2]=1;
00401 triggconst[6][3]=0.426516;
00402 triggconst[6][4]=0.828706;
00403 triggconst[6][5]=1;
00404 triggconst[6][6]=0.865995;
00405 triggconst[6][7]=1; 
00406 
00407 triggconst[7][0]=1;
00408 triggconst[7][1]=1;
00409 triggconst[7][2]=1;
00410 triggconst[7][3]=1;
00411 triggconst[7][4]=1;
00412 triggconst[7][5]=1;
00413 triggconst[7][6]=1;
00414 triggconst[7][7]=1; 
00415 
00416 }
00417 
00418 
00419 
00420 Float_t TTriggerReweight::trigeff_9899__(Float_t *x, Float_t *y, Int_t *imeth, Int_t *iyear)
00421 {
00422 
00423     static Float_t pa[9] = { (Float_t)1.,(Float_t)1.,(Float_t)1.,(Float_t)1.,(Float_t)1.,(
00424             float)1.,(Float_t)1.,(Float_t)1.,(Float_t)1. };
00425     static Float_t pc4[9] = { (Float_t)8.48,(Float_t)8.51,(Float_t)4.89,(Float_t)4.29,(
00426             float)3.9,(Float_t)6.16,(Float_t)7.84,(Float_t)10.91,(Float_t)11.94 };
00427     static Float_t pb5[9] = { (Float_t)36.55,(Float_t)9.22,(Float_t)8.54,(Float_t)9.99,(
00428             float)9.77,(Float_t)9.01,(Float_t)8.82,(Float_t)9.42,(Float_t)9.52 };
00429     static Float_t pc5[9] = { (Float_t)8.97,(Float_t)6.07,(Float_t)4.93,(Float_t)2.73,(
00430             float)3.46,(Float_t)5.02,(Float_t)7.29,(Float_t)8.38,(Float_t)8.39 };
00431     static Float_t pb6[9] = { (Float_t)36.55,(Float_t)9.17,(Float_t)8.29,(Float_t)9.8,(
00432             float)9.7,(Float_t)8.78,(Float_t)8.9,(Float_t)9.31,(Float_t)8.92 };
00433     static Float_t pc6[9] = { (Float_t)8.97,(Float_t)6.07,(Float_t)4.99,(Float_t)2.64,(
00434             float)3.13,(Float_t)5.01,(Float_t)7.03,(Float_t)8.28,(Float_t)8.17 };
00435     static Float_t pb7[9] = { (Float_t)36.55,(Float_t)9.22,(Float_t)8.6,(Float_t)9.6,(
00436             float)9.24,(Float_t)7.97,(Float_t)7.88,(Float_t)9.04,(Float_t)8.62 };
00437     static Float_t pc7[9] = { (Float_t)8.97,(Float_t)6.07,(Float_t)4.85,(Float_t)2.76,(
00438             float)3.26,(Float_t)5.28,(Float_t)7.48,(Float_t)8.04,(Float_t)7.3 };
00439     static Float_t pb8[9] = { (Float_t)36.55,(Float_t)9.17,(Float_t)8.34,(Float_t)9.45,(
00440             float)9.34,(Float_t)8.,(Float_t)7.98,(Float_t)9.03,(Float_t)8.19 };
00441     static Float_t pc8[9] = { (Float_t)8.97,(Float_t)6.07,(Float_t)4.91,(Float_t)2.68,(
00442             float)2.99,(Float_t)5.15,(Float_t)7.29,(Float_t)7.91,(Float_t)7.31 };
00443     static Float_t thm1[9] = { (Float_t)9.34,(Float_t)12.64,(Float_t)17.36,(Float_t)24.8,
00444             (Float_t)34.73,(Float_t)48.79,(Float_t)69.27,(Float_t)93.48,(Float_t)120.86 
00445             };
00446     static Float_t thm2[9] = { (Float_t)9.32,(Float_t)12.63,(Float_t)17.39,(Float_t)
00447             24.84,(Float_t)34.78,(Float_t)48.77,(Float_t)69.24,(Float_t)93.68,(Float_t)
00448             122.5 };
00449     static Float_t thm3[9] = { (Float_t)9.34,(Float_t)12.64,(Float_t)17.37,(Float_t)
00450             24.96,(Float_t)34.79,(Float_t)48.79,(Float_t)69.31,(Float_t)93.73,(Float_t)
00451             122.15 };
00452     static Float_t thm4[9] = { (Float_t)9.32,(Float_t)12.63,(Float_t)17.4,(Float_t)24.92,
00453             (Float_t)34.8,(Float_t)48.8,(Float_t)69.31,(Float_t)93.84,(Float_t)122.58 };
00454     static Float_t thm5[9] = { (Float_t)9.33,(Float_t)12.63,(Float_t)17.4,(Float_t)24.83,
00455             (Float_t)34.76,(Float_t)48.83,(Float_t)69.64,(Float_t)93.44,(Float_t)120.93 
00456             };
00457     static Float_t thm6[9] = { (Float_t)9.33,(Float_t)12.63,(Float_t)17.41,(Float_t)
00458             24.84,(Float_t)34.76,(Float_t)48.83,(Float_t)69.59,(Float_t)93.44,(Float_t)
00459             121.57 };
00460     static Float_t thm7[9] = { (Float_t)9.33,(Float_t)12.63,(Float_t)17.41,(Float_t)
00461             24.86,(Float_t)34.77,(Float_t)48.88,(Float_t)69.59,(Float_t)93.49,(Float_t)
00462             121.62 };
00463     static Float_t thm8[9] = { (Float_t)9.33,(Float_t)12.63,(Float_t)17.41,(Float_t)
00464             24.86,(Float_t)34.77,(Float_t)48.87,(Float_t)69.58,(Float_t)93.52,(Float_t)
00465             121.87 };
00466     static Float_t pb1[9] = { (Float_t)39.32,(Float_t)10.01,(Float_t)10.72,(Float_t)
00467             10.61,(Float_t)10.53,(Float_t)10.35,(Float_t)10.2,(Float_t)9.96,(Float_t)
00468             10.4 };
00469     static Float_t pc1[9] = { (Float_t)8.36,(Float_t)8.53,(Float_t)4.97,(Float_t)4.86,(
00470             float)4.9,(Float_t)6.82,(Float_t)8.79,(Float_t)12.93,(Float_t)15.19 };
00471     static Float_t pb2[9] = { (Float_t)38.92,(Float_t)9.84,(Float_t)10.21,(Float_t)9.58,(
00472             float)9.85,(Float_t)9.41,(Float_t)9.86,(Float_t)9.4,(Float_t)8.5 };
00473     static Float_t pc2[9] = { (Float_t)8.48,(Float_t)8.51,(Float_t)4.92,(Float_t)4.58,(
00474             float)4.19,(Float_t)6.5,(Float_t)8.11,(Float_t)11.95,(Float_t)12.74 };
00475     static Float_t pb3[9] = { (Float_t)39.32,(Float_t)10.01,(Float_t)10.7,(Float_t)10.01,
00476             (Float_t)9.5,(Float_t)8.75,(Float_t)9.32,(Float_t)9.03,(Float_t)7.32 };
00477     static Float_t pc3[9] = { (Float_t)8.36,(Float_t)8.53,(Float_t)4.88,(Float_t)4.5,(
00478             float)4.18,(Float_t)6.49,(Float_t)8.04,(Float_t)10.99,(Float_t)12.52 };
00479     static Float_t pb4[9] = { (Float_t)38.92,(Float_t)9.84,(Float_t)10.19,(Float_t)9.17,(
00480             float)9.09,(Float_t)8.57,(Float_t)9.16,(Float_t)8.61,(Float_t)6.13 };
00481 
00482     /* System generated locals */
00483     Float_t ret_val;
00484 
00485     /* Local variables */
00486     static Int_t i__;
00487     static Float_t pb[9], pc[9], ef1, ef2, thm[9];
00488 
00489 
00490 /* x (pt_had_AEFR) */
00491 /* y (theta_had_AEFR) */
00492 /* imeth=1(st 66 77), =2(st 66 67 77), =3(st 66 71 77), =4(st 66 67 71 77) */
00493 /* iyear=98(Etmiss t0-1), 99(Etmiss t0) */
00494 /*      dimension pa(9),pb(9),pc(9),thm(9) */
00495     ret_val = (Float_t)1.;
00496     if (*iyear != 98 && *iyear != 99) {
00497         printf("check the input year (98 or 99)\n");
00498         return ret_val;
00499     }
00500     if (*iyear == 98) {
00501         for (i__ = 1; i__ <= 9; ++i__) {
00502             if (*imeth == 1) {
00503                 pb[i__ - 1] = pb1[i__ - 1];
00504                 pc[i__ - 1] = pc1[i__ - 1];
00505                 thm[i__ - 1] = thm1[i__ - 1];
00506             } else if (*imeth == 2) {
00507                 pb[i__ - 1] = pb2[i__ - 1];
00508                 pc[i__ - 1] = pc2[i__ - 1];
00509                 thm[i__ - 1] = thm2[i__ - 1];
00510             } else if (*imeth == 3) {
00511                 pb[i__ - 1] = pb3[i__ - 1];
00512                 pc[i__ - 1] = pc3[i__ - 1];
00513                 thm[i__ - 1] = thm3[i__ - 1];
00514             } else if (*imeth == 4) {
00515                 pb[i__ - 1] = pb4[i__ - 1];
00516                 pc[i__ - 1] = pc4[i__ - 1];
00517                 thm[i__ - 1] = thm4[i__ - 1];
00518             }
00519         }
00520     } else if (*iyear == 99) {
00521         for (i__ = 1; i__ <= 9; ++i__) {
00522             if (*imeth == 1) {
00523                 pb[i__ - 1] = pb5[i__ - 1];
00524                 pc[i__ - 1] = pc5[i__ - 1];
00525                 thm[i__ - 1] = thm5[i__ - 1];
00526             } else if (*imeth == 2) {
00527                 pb[i__ - 1] = pb6[i__ - 1];
00528                 pc[i__ - 1] = pc6[i__ - 1];
00529                 thm[i__ - 1] = thm6[i__ - 1];
00530             } else if (*imeth == 3) {
00531                 pb[i__ - 1] = pb7[i__ - 1];
00532                 pc[i__ - 1] = pc7[i__ - 1];
00533                 thm[i__ - 1] = thm7[i__ - 1];
00534             } else if (*imeth == 4) {
00535                 pb[i__ - 1] = pb8[i__ - 1];
00536                 pc[i__ - 1] = pc8[i__ - 1];
00537                 thm[i__ - 1] = thm8[i__ - 1];
00538             }
00539         }
00540     }
00541     if (*y < thm[0]) {
00542         if (*y < (Float_t)10.) {
00543             ef1 = pa[0] / (exp((pb[0] - *x) / pc[0]) + (Float_t)1.);
00544         } else {
00545             ef1 = pa[0] * ((Float_t)1. - exp((pb[0] - *x) / pc[0]));
00546         }
00547         ret_val = ef1;
00548     } else if (*y > thm[8]) {
00549         ef2 = pa[8] * ((Float_t)1. - exp((pb[8] - *x) / pc[8]));
00550         ret_val = ef2;
00551     } else {
00552         for (i__ = 1; i__ <= 8; ++i__) {
00553             if (*y > thm[i__ - 1] && *y < thm[i__]) {
00554                 if (i__ == 1) {
00555                     ef1 = pa[i__ - 1] / (exp((pb[i__ - 1] - *x) / pc[i__ - 1])
00556                              + (Float_t)1.);
00557                 } else {
00558                     ef1 = pa[i__ - 1] * ((Float_t)1. - exp((pb[i__ - 1] - *x) / 
00559                             pc[i__ - 1]));
00560                 }
00561                 ef2 = pa[i__] * ((Float_t)1. - exp((pb[i__] - *x) / pc[i__]));
00562                 ret_val = ef1 + (*y - thm[i__ - 1]) / (thm[i__] - thm[i__ - 1]
00563                         ) * (ef2 - ef1);
00564                 goto L999;
00565             }
00566         }
00567     }
00568 L999:
00569     if (ret_val < (Float_t)0.) {
00570         ret_val = (Float_t)0.;
00571     }
00572     if (ret_val > (Float_t)1.) {
00573         ret_val = (Float_t)1.;
00574     }
00575     return ret_val;
00576 }
00577 
00578 
00579 Float_t TTriggerReweight::trig_corr__(Float_t *pet, Float_t *gar, Int_t *iyr)
00580 {
00581 
00582 static Int_t c__1 = 1;
00583 static Int_t c__2 = 2;
00584 static Int_t c__3 = 3;
00585 static Int_t c__4 = 4;
00586 static Int_t c__5 = 5;
00587 static Int_t c__6 = 6;
00588 static Int_t c__7 = 7;
00589 static Int_t c__8 = 8;
00590 
00591 
00592     /* System generated locals */
00593     Float_t ret_val = 0;
00594 
00595     if (*iyr == 99) {
00596       //      printf("*gar = %f\n",*gar);
00597       //      printf("gbins0 = %f gbins1 = %f gbins2 = %f gbins3 = %f gbins4 = %f gbins5 = %f gbins6 = %f gbins7 = %f gbins8 = %f\n",tr99_1.gbins[0],tr99_1.gbins[1],tr99_1.gbins[2],tr99_1.gbins[3],tr99_1.gbins[4],tr99_1.gbins[5],tr99_1.gbins[6],tr99_1.gbins[7],tr99_1.gbins[8]);
00598         if (*gar >= tr99_1.gbins[0] && *gar < tr99_1.gbins[1]) {
00599             fill_(&c__1, iyr);
00600             ret_val = tr99_1.a[0] - tr99_1.b[0] * exp(-tr99_1.c__[0] * *pet);
00601         }
00602         if (*gar >= tr99_1.gbins[1] && *gar < tr99_1.gbins[2]) {
00603             inter_(&c__2, &c__3, gar, pet, &c__2, iyr);
00604             ret_val = tr99_1.tr_val__[0];
00605         }
00606         if (*gar >= tr99_1.gbins[2] && *gar < tr99_1.gbins[3]) {
00607             inter_(&c__3, &c__4, gar, pet, &c__3, iyr);
00608             ret_val = tr99_1.tr_val__[0];
00609         }
00610         if (*gar >= tr99_1.gbins[3] && *gar < tr99_1.gbins[4]) {
00611             inter_(&c__4, &c__5, gar, pet, &c__4, iyr);
00612             ret_val = tr99_1.tr_val__[0];
00613         }
00614         if (*gar >= tr99_1.gbins[4] && *gar < tr99_1.gbins[5]) {
00615             inter_(&c__5, &c__6, gar, pet, &c__5, iyr);
00616             ret_val = tr99_1.tr_val__[0];
00617         }
00618         if (*gar >= tr99_1.gbins[5] && *gar < tr99_1.gbins[6]) {
00619             inter_(&c__6, &c__7, gar, pet, &c__6, iyr);
00620             ret_val = tr99_1.tr_val__[0];
00621         }
00622         if (*gar >= tr99_1.gbins[6] && *gar < tr99_1.gbins[7]) {
00623             inter_(&c__7, &c__8, gar, pet, &c__7, iyr);
00624             ret_val = tr99_1.tr_val__[0];
00625         }
00626         if (*gar >= tr99_1.gbins[7] && *gar < tr99_1.gbins[8]) {
00627             fill_(&c__8, iyr);
00628             ret_val = tr99_1.a[0] - tr99_1.b[0] * exp(-tr99_1.c__[0] * *pet);
00629         }
00630     }
00631     if (ret_val > (float)1.) {
00632         ret_val = (float)1.;
00633     }
00634     if (ret_val < (float)0.) {
00635         ret_val = (float)0.;
00636     }
00637     return ret_val;
00638 } /* trig_corr__ */
00639 
00640 
00641 Int_t TTriggerReweight::fill_(Int_t *sel, Int_t *iyr)
00642 {
00643 // Need for trigger efficienies calculation
00644 // Fills common block 
00645     if (*iyr == 99) {
00646         tr99_1.a[0] = tr99_1.out3[*sel * 3 - 3];
00647         tr99_1.b[0] = tr99_1.out3[*sel * 3 - 2];
00648         tr99_1.c__[0] = tr99_1.out3[*sel * 3 - 1];
00649     }
00650     return 0;
00651 
00652 }
00653 
00654 Int_t TTriggerReweight::inter_(Int_t *g1,Int_t *g2,Float_t *gar,Float_t *pet,Int_t *sel, Int_t *iyr){
00655 // For trigger efficiencies
00656 // Interpolate between 2 bins
00657 
00658     /* Local variables */
00659     static Float_t ratio, lower, upper, al, am, bm, cm, bl, cl, au, bu, cu,
00660             middle, splitl, splitm, splitu, inv;
00661 
00662 /*     Routine to interpolate between gamma-binned efficiencies */
00663 /*     Particular gamma value gar lies between g1 and g2. */
00664     if (*iyr == 99) {
00665         am = tr99_1.out3[*sel * 3 - 3];
00666         bm = tr99_1.out3[*sel * 3 - 2];
00667         cm = tr99_1.out3[*sel * 3 - 1];
00668         al = tr99_1.out3[(*sel - 1) * 3 - 3];
00669         bl = tr99_1.out3[(*sel - 1) * 3 - 2];
00670         cl = tr99_1.out3[(*sel - 1) * 3 - 1];
00671         au = tr99_1.out3[(*sel + 1) * 3 - 3];
00672         bu = tr99_1.out3[(*sel + 1) * 3 - 2];
00673         cu = tr99_1.out3[(*sel + 1) * 3 - 1];
00674         lower = al - bl * exp(-cl * *pet);
00675         middle = am - bm * exp(-cm * *pet);
00676         upper = au - bu * exp(-cu * *pet);
00677         splitl = (tr99_1.gbins[*g1 - 2] + tr99_1.gbins[*g1 - 1]) / 2;
00678         splitm = (tr99_1.gbins[*g2 - 1] + tr99_1.gbins[*g1 - 1]) / 2;
00679         splitu = (tr99_1.gbins[*g2 - 1] + tr99_1.gbins[*g2]) / 2;
00680     }
00681     if (*gar >= splitm) {
00682         ratio = (*gar - splitm) / (splitu - splitm);
00683         inv = 1 - ratio;
00684         tr99_1.tr_val__[0] = ratio * upper + inv * middle;
00685     }
00686     if (*gar < splitm) {
00687         ratio = (splitm - *gar) / (splitm - splitl);
00688         inv = 1 - ratio;
00689         tr99_1.tr_val__[0] = ratio * lower + inv * middle;
00690     }
00691     return 0;
00692 
00693 }