Attachment 'CreateDataForSignalAndBackground.C'

Download

   1   #include "WIMP_functions.C"
   2 
   3 #include <fstream>
   4 #include "TFile.h"
   5 #include "TProfile.h"
   6 
   7 
   8 
   9 //extern "C" void differentialxsection_(double* s, double* X_photon, double* Theta_photon, double* P_e, double* P_p, double* sigma_diff);
  10 
  11 
  12 
  13 void CreateDataForSignalAndBackground(){
  14   //void CreateDataForSignalAndBackground(double l, double pe, double pp){
  15 
  16   // ####### PARAMETERS
  17   
  18   //double lumi = l;
  19   //double Pe = pe;
  20   //double Pp = pp;
  21   
  22   double lumi = 2000.0;
  23   double Pe = 0.0;
  24   double Pp = 0.0;
  25   
  26   double cos_cut = 0.995; //cos_cut = 0.98
  27       
  28   //_lambda=(double)100;// p-wave = 15, s-wave = 32
  29 
  30   TFile* f;
  31   char filename[256];
  32   char name[256];
  33   char obj_title[256];
  34 
  35 
  36   // #########################################################################
  37   // #########################################################################
  38   // ##################   BACKGROUND DATA TEMPLATES ######################
  39   // #########################################################################
  40   
  41   //##################### OUTPUTFILES #####################
  42 //  TFile* f_data_bg[7];
  43 //  TH1D* h_data_bg_E[7];
  44 //  
  45 //  for (int m(0); m < 7 ; m++){
  46 //    if (m == 0) {
  47 //      Pe = 0.0;
  48 //      Pp = 0.0;
  49 //    }
  50 //    if (m == 1) {
  51 //      Pe = 0.8;
  52 //      Pp = 0.0;
  53 //    }
  54 //    if (m == 2) {
  55 //      Pe = 0.8;
  56 //      Pp = -0.3;
  57 //    }
  58 //    if (m == 3) {
  59 //      Pe = 0.8;
  60 //      Pp = -0.6;
  61 //    }
  62 //    if (m == 4) {
  63 //      Pe = -0.8;
  64 //      Pp = 0.0;
  65 //    }
  66 //    if (m == 5) {
  67 //      Pe = -0.8;
  68 //      Pp = 0.3;
  69 //    }
  70 //    if (m == 6) {
  71 //      Pe = -0.8;
  72 //      Pp = 0.6;
  73 //    }
  74 //  
  75 //   
  76 //    sprintf(filename,"../data/Templates/Data_Background_Pe%.1f_Pp%.1f_Lumi%.0f.root",Pe,Pp,lumi);
  77 //    f_data_bg[m] = new TFile(filename,"RECREATE");
  78 //      
  79 //    sprintf(obj_title,"Background_Pe%.1f_Pp%.1f_Lumi%.0f",Pe,Pp,lumi);
  80 //    h_data_bg_E[m] = new TH1D(obj_title,"Background;E [GeV]; Entries",300,0.0,300.0);
  81 //    h_data_bg_E[m]->Sumw2();
  82 //  }
  83 //   
  84 //  
  85 //  // ################ CREATE BACKGROUND TEMPLATES ############
  86 //    
  87 //  ifstream txtfile("../steeringfiles/ILD00_Data_RootFiles.txt");
  88 //   
  89 //    
  90 //  while (!txtfile.eof()){
  91 //    txtfile.getline(name,256,'\n');
  92 //    printf("%s\n",name);
  93 //    f = new TFile(name,"OPEN");
  94 //    if (f->IsZombie()) continue;
  95 //      
  96 //    // ********** GET BACKGROUND DATA TREES
  97 //      
  98 //    TTree* tree = (TTree*)f->Get("Tree");
  99 //    TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 100 //    TTree* globaltree = (TTree*)f->Get("GlobalParameters");
 101 //    LoadTree(tree);
 102 //    LoadSelPhotonTree(seltree);
 103 //    LoadGlobalTree(globaltree);
 104 //    
 105 //    globaltree->GetEntry(0);
 106 //    
 107 //    double factor = (double)NumberOfEventsInFile/NumberOfEventsInRun;
 108 //    // ********** LOOP OVER ENTRIES
 109 //       
 110 //    int n_entries = seltree->GetEntries();
 111 //    //n_entries = (int)n_entries/1000;
 112 //      
 113 //    for (int j(0); j < n_entries; j++){
 114 //      if (j%10000 == 0) printf("event %i\n",j);
 115 //      tree->GetEntry(j);
 116 //      seltree->GetEntry(j);
 117 //      MCLumi_Run = MCLumi_Run*factor;
 118 //      for (int m(0); m < 7 ; m++){
 119 //  	if (m == 0) {
 120 //  	  Pe = 0.0;
 121 //  	  Pp = 0.0;
 122 //  	}
 123 //  	if (m == 1) {
 124 //  	  Pe = 0.8;
 125 //  	  Pp = 0.0;
 126 //  	}
 127 //  	if (m == 2) {
 128 //  	  Pe = 0.8;
 129 //  	  Pp = -0.3;
 130 //  	}
 131 //  	if (m == 3) {
 132 //  	  Pe = 0.8;
 133 //  	  Pp = -0.6;
 134 //  	}
 135 //  	if (m == 4) {
 136 //  	  Pe = -0.8;
 137 //  	  Pp = 0.0;
 138 //  	}
 139 //  	if (m == 5) {
 140 //  	  Pe = -0.8;
 141 //  	  Pp = 0.3;
 142 //  	}
 143 //  	if (m == 6) {
 144 //  	  Pe = -0.8;
 145 //  	  Pp = 0.6;
 146 //  	}
 147 //  	CreatePolarisationWeights(Pe,Pp);
 148 //  	double weight = GetEventWeight(lumi, P_e, P_p)*EventWeight;
 149 //  	h_data_bg_E[m]->Fill(SelPhoton_energy[0],weight);
 150 //      }
 151 //    }
 152 //    delete f;
 153 //  }
 154 //  
 155 //  // #########  ADD BHABHA EXPECTATION HISTOGRAM
 156 //  char polstate[128];
 157 //  char histoname[128];
 158 //  TFile* f_bhabha = new TFile("../data/Templates/BhabhaBackground.root","OPEN");
 159 //  
 160 //  for ( int m(0); m < 7; m++){
 161 //    switch(m){
 162 //    case 0:
 163 //      sprintf(polstate,"p0p0");
 164 //      break;
 165 //    case 1:
 166 //      sprintf(polstate,"p8p0");
 167 //      break;
 168 //    case 2:
 169 //      sprintf(polstate,"p8m3");
 170 //      break;
 171 //    case 3:
 172 //      sprintf(polstate,"p8m6");
 173 //      break;
 174 //    case 4:
 175 //      sprintf(polstate,"m8p0");
 176 //      break;
 177 //    case 5:
 178 //      sprintf(polstate,"m8p3");
 179 //      break;
 180 //    case 6:
 181 //      sprintf(polstate,"m8p6");
 182 //      break;
 183 //    }
 184 //    sprintf(histoname,"Bhabha_E_Lumi50_%s",polstate);
 185 //    TH1D* h_bhabha_E = (TH1D*)f_bhabha->Get(histoname);
 186 //    h_bhabha_E->Scale(lumi/50.0);
 187 //    h_data_bg_E[m]->Add(h_bhabha_E,1.0);
 188 //  }
 189 //  
 190 //  for ( int m(0); m < 7; m++){
 191 //    f_data_bg[m]->cd();
 192 //    h_data_bg_E[m]->Write();
 193 //    f_data_bg[m]->Close();
 194 //  }
 195 //  f_bhabha->Close();
 196 //    
 197 //  // ###### FOR COMPARISON, STORE PARAMETRISED BACKGROUND
 198 // // TFile* f_full_bg = new TFile("../data/Templates/ParametrisedBackground.root","OPEN");
 199 // // sprintf(polstate,"m830");  
 200 // // sprintf(obj_title,"Background_E_Modeled_Lumi50_%s",polstate);
 201 // // TH1D* h_full_bg_template = (TH1D*)f_full_bg->Get(obj_title);
 202 // // h_full_bg_template->Scale(lumi/50.0);
 203 //
 204 //
 205 //
 206 //
 207 //  
 208 
 209 
 210   // #########################################################################
 211   // #########################################################################
 212   // ##################   SIGNAL DATA TEMPLATES ############################
 213   // #########################################################################
 214   
 215 
 216 
 217   FILE* wimpfile = fopen("../steeringfiles/WIMP.steer","r");
 218   LoadWIMPParameters(wimpfile);
 219   _BR = 1.0;
 220   _BR_eff = 1.0;
 221   
 222   int _sigma_eRpL = 0.0;
 223   int _sigma_eRpR = 0.0;
 224   int _sigma_eLpL = 0.0;
 225   int _sigma_eLpR = 0.0;
 226   
 227 //   int k= 0; // Vector	    
 228   int k= 1; // Scalar,s-channel
 229 //  int k= 2; // Axial-vector 
 230   
 231   //##################### OUTPUTFILES #####################
 232 
 233 
 234   char _operator[20];
 235  
 236   TFile* f_sig[9];
 237   TProfile* h_weight[8];
 238   TH1D* h_array[9][250][10][1]; // pol, M , lambda ,k-operator
 239  
 240   for (int m(0); m < 9 ; m++){
 241 
 242     if (m == 0) {
 243       Pe = 0.0;
 244       Pp = 0.0;
 245     }
 246     if (m == 1) {
 247       Pe = 0.8;
 248       Pp = 0.0;
 249     }
 250     if (m == 2) {
 251       Pe = 0.8;
 252       Pp = 0.3;
 253     }
 254     if (m == 3) {
 255       Pe = 0.8;
 256       Pp = 0.6;
 257     }
 258     if (m == 4) {
 259       Pe = -0.8;
 260       Pp = 0.0;
 261     } 	
 262     if (m == 5) {
 263       Pe = -0.8;
 264       Pp = 0.3;
 265     }
 266     if (m == 6) {
 267       Pe = -0.8;
 268       Pp = 0.6;
 269     }
 270     if (m == 7) {
 271       Pe = 0.8;
 272       Pp = -0.3;
 273     }
 274     if (m == 8) {
 275       Pe = 0.8;
 276       Pp = -0.6;
 277     }
 278     
 279     sprintf(filename,"../data/Templates_Andrii/temp/1D_Data_Signal_Pe%.1f_Pp%.1f_Lumi%.0f.root",Pe,Pp,lumi);
 280     f_sig[m] = new TFile(filename,"RECREATE");
 281    
 282     h_weight[m] = new TProfile("weights","weights;cos #Theta;w",100,-1.0,1.0);
 283 
 284     for (int i(0); i < 250; i++){
 285       for (int l(0); l< 1; l++){//l = Lambda
 286 
 287 	
 288 	    if (k == 0){
 289 	      sprintf(_operator,"Vector");
 290 	    }
 291 	    if (k == 1){
 292 	      sprintf(_operator,"Scalar,s-channel");
 293 	    }
 294 	    if (k == 2){
 295 		sprintf(_operator,"Axial-vector");
 296 	    }	   
 297 	    if (k == 3){
 298 		sprintf(_operator,"Scalar,t-channel");
 299 	    }
 300 	    if(i==0)
 301 	    cout<<_operator<<endl;	    
 302 	    
 303 	    _mass = (double)i + 1.0; //100.0;
 304 	    _spin= 1.0/2.0;
 305 	    _lambda=(double)(l+1)*1580;	    
 306 
 307 	    sprintf(obj_title,"M%.0f_S%.1f_%s_Lambda%.0f",_mass,_spin,_operator,_lambda);
 308 	    h_array[m][i][l][0] = new TH1D(obj_title,obj_title,300,0.0,300.0);
 309 	    h_array[m][i][l][0]->Sumw2();
 310 	  
 311       }
 312     }
 313   }
 314 
 315 
 316   double weight;
 317   double xsec_bg;
 318   double xsec_sig;
 319 
 320   // ---------------------------------------------------------
 321 
 322   // ###################################################################
 323   // ################ CREATE SIGNAL TEMPLATES ##########################
 324   // ################### METHOD 1, DIRECT CALCULATION OF WEIGHTS #######
 325  
 326   //  CreatePolarisationWeights(0.0,0.0);
 327   // 
 328   //
 329   //  //----------------------------------------------------------------------------  
 330   //
 331   //  ifstream txtfile2("../steeringfiles/ILD00_Signal_RootFiles.txt");
 332   //
 333   //  while (!txtfile2.eof()){
 334   //    txtfile2.getline(name,256,'\n');
 335   //   
 336   //    printf("%s\n",name);
 337   //    f = new TFile(name,"OPEN");
 338   //    if (f->IsZombie()) continue;
 339   //   
 340   //    // ********** GET BACKGROUND DATA TREES
 341   //   
 342   //    TTree* tree = (TTree*)f->Get("Tree");
 343   //    TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 344   //    TTree* xsectree = (TTree*)f->Get("XSecTree");
 345   //    LoadTree(tree);
 346   //    LoadSelPhotonTree(seltree);
 347   //    LoadXSecTree(xsectree);
 348   //   
 349   //    // #########  LOOP OVER ENTRIES
 350   //    int nentries = seltree->GetEntries();
 351   //    // nentries = 1e2;
 352   //
 353   //    for (int l(0); l < nentries; l++){
 354   //      if (l%10000 == 0) printf("event %i\n",l);      
 355   //
 356   //      tree->GetEntry(l);
 357   //      seltree->GetEntry(l);
 358   //      xsectree->GetEntry(l);
 359   //           
 360   //      weight = GetEventWeight(lumi,MCLumi_Run,P_e,P_p);
 361   //      double x  = E_prime*2.0/sqrt(1D_Data_Signal_Pe0.0_Pp0.0_Lumi50.roots_prime);
 362   //      Pe = 0.0;
 363   //      Pp = 0.0;
 364   //      differentialxsection_(&s_prime,&x, &Theta_prime,&Pe,&Pp,&xsec_bg); //dsigma/dx/dcos(theta) in pb 
 365   //      // s_prime = 500.0*500.0;
 366   //      // x = SelPhoton_energy[0]*2.0/sqrt(s_prime);
 367   //      // Theta_prime = TMath::ACos(SelPhoton_cosTheta[0]);
 368   //      // differentialxsection_(&s_prime,&x, &Theta_prime,&Pe,&Pp,&xsec_bg);
 369   //
 370   //
 371   //      for (int m(0); m < 7 ; m++){
 372   //	if (m == 0) {
 373   //	  Pe = 0.0;
 374   //	  Pp = 0.0;
 375   //	}
 376   //	if (m == 1) {
 377   //	  Pe = 0.8;
 378   //	  Pp = 0.0;
 379   //	}
 380   //	if (m == 2) {
 381   //	  Pe = 0.8;
 382   //	  Pp = -0.3;
 383   //	}
 384   //	if (m == 3) {
 385   //	  Pe = 0.8;
 386   //	  Pp = -0.6;
 387   //	}
 388   //	if (m == 4) {
 389   //	  Pe = -0.8;
 390   //	  Pp = 0.0;
 391   //	}
 392   //	if (m == 5) {
 393   //	  Pe = -0.8;
 394   //	  Pp = 0.3;
 395   //	}
 396   //	if (m == 6) {
 397   //	  Pe = -0.8;
 398   //	  Pp = 0.6;
 399   //	}
 400   //
 401   //	for (int i(0); i < 150; i++){
 402   //	  for (int s(0); s < 3; s++){
 403   //	    for ( int j(0); j < 2; j++){
 404   //	      for ( int k(0); k < 3; k++){
 405   //		if (k == 0){
 406   //		  _kappa_eRpL = 1.0;
 407   //		  _kappa_eRpR = 1.0;
 408   //		  _kappa_eLpL = 1.0;
 409   //		  _kappa_eLpR = 1.0;
 410   //		}
 411   //		if (k == 1){
 412   //		  _kappa_eRpL = 2.0;
 413   //		  _kappa_eRpR = 0.0;
 414   //		  _kappa_eLpL = 0.0;
 415   //		  _kappa_eLpR = 2.0;
 416   //		}
 417   //		if (k == 2){
 418   //		  _kappa_eRpL = 4.0;
 419   //		  _kappa_eRpR = 0.0;
 420   //		  _kappa_eLpL = 0.0;
 421   //		  _kappa_eLpR = 0.0;
 422   //		}
 423   //		_mass = (double)i + 100.0;
 424   //		_spin = (double)s/2.0;
 425   //		_annihilator = (double)j;
 426   //		if (j == 1) _sigma_an = 6.0;
 427   //		else   _sigma_an = 0.85; 
 428   //
 429   //		_BR_eff = _BR * GetWIMP_PolFactor(_kappa_eRpL, _kappa_eRpR, _kappa_eLpL,_kappa_eLpR,Pe,Pp); 
 430   //		xsec_sig = DiffWIMPCrossSection(sqrt(s_prime), E_prime, Theta_prime);  //dsigma/dx/dcos(theta) in pb 
 431   //		// xsec_sig = DiffWIMPCrossSection(sqrt(s_prime),SelPhoton_energy[0] , Theta_prime);
 432   //
 433   //		if (xsec_bg > 0.0){
 434   //		  double weight2 = weight * xsec_sig/xsec_bg;
 435   //
 436   //		  h_weight[m]->Fill(TMath::Cos(Theta_prime),weight2);
 437   //		  h_array[m][i][s][j][k]->Fill(SelPhoton_energy[0], weight2);
 438   //		}
 439   //	      }
 440   //	    }
 441   //	  }
 442   //	}
 443   //      }
 444   //    }
 445   //    delete f; 
 446   //  }
 447 
 448 
 449 
 450 
 451   // ###################################################################
 452   // ################ CREATE SIGNAL TEMPLATES ##########################
 453   // ################### METHOD 2, Weights with dsig/de #######
 454 
 455 
 456   // ########### OPEN BACKGROUND CROSS SECTION FILE
 457   sprintf(filename,"../data/Auxiliary/BackgroundExpectationPe0.0Pp0.0.root");
 458   TFile* f_bg_exp = new TFile(filename,"OPEN");
 459   TH2D* h_bg_2d = (TH2D*)f_bg_exp->Get("nunugamma2d"); 
 460   
 461   
 462   // ######## PROJECTION ON ENERGY 
 463     
 464   double w = h_bg_2d->GetYaxis()->GetBinWidth(1);
 465   int binx = h_bg_2d->GetXaxis()->GetNbins();
 466 
 467   int biny = h_bg_2d->GetYaxis()->GetNbins();
 468 
 469   TH1D* h_bg_dsigmadE = new TH1D("nunugammaE","nunugammaE",binx,0.0,300.0); // dsigma/dE in fb
 470     
 471   for ( int i(1); i < binx+1; i++){
 472     double cont(0.0);
 473     for ( int j(1); j < biny+1; j++){
 474        if ( TMath::Abs(h_bg_2d->GetYaxis()->GetBinCenter(j)) < cos_cut){
 475   	cont += h_bg_2d->GetBinContent(i,j) * w;
 476       }
 477       h_bg_dsigmadE->SetBinContent(i,cont);
 478     }
 479   }
 480 
 481 
 482   _BR_eff = 1.0;
 483   TH1D* h_WIMPdsigmadE[250][1][10]; // dsigma/dE in fb
 484   for (int i(0); i < 250; i++){
 485 
 486       int s(0);
 487 
 488       	for (int l(0);l<1;l++){
 489 	
 490 	_mass = (double)i + 1.0; // 100.0;
 491   	_spin = (double)1/2.0;
 492 	_lambda=(double)(l+1)*1580;
 493  
 494   	h_WIMPdsigmadE[i][s][l]= new TH1D(Form("DsigmadE_WIMP%i%i%i",i,s,l),"DsigmadE_WIMP;E[GeV];d#sigma /dE [fb]",300,0.0,300.0);
 495   	DsigmadE_WIMP_new(h_WIMPdsigmadE[i][s][l],500.0,cos_cut);
 496 
 497 	}
 498 
 499   }
 500       
 501     
 502   CreatePolarisationWeights(0.0,0.0);
 503    
 504   
 505   //----------------------------------------------------------------------------  
 506   
 507   ifstream txtfile2("../steeringfiles/ILD00_Signal_RootFiles.txt");
 508   
 509   while (!txtfile2.eof()){
 510     txtfile2.getline(name,256,'\n');
 511      
 512     printf("%s\n",name);
 513     f = new TFile(name,"OPEN");
 514     if (f->IsZombie()) continue;
 515      
 516     // ********** GET BACKGROUND DATA TREES
 517      
 518     TTree* tree = (TTree*)f->Get("Tree");
 519     TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 520     TTree* globaltree = (TTree*)f->Get("GlobalParameters");
 521      // TTree* xsectree = (TTree*)f->Get("XSecTree");
 522     LoadTree(tree);
 523     LoadSelPhotonTree(seltree);
 524     // LoadXSecTree(xsectree);
 525     LoadGlobalTree(globaltree);
 526     globaltree->GetEntry(0);
 527 
 528     double factor = (double)NumberOfEventsInFile/NumberOfEventsInRun;  
 529     // #########  LOOP OVER ENTRIES
 530     int nentries = seltree->GetEntries();
 531 //       nentries = (int)1e0;
 532   
 533     for (int l(0); l < nentries; l++){
 534       if (l%1000 == 0){ printf("event %i\n",l);fflush(stdout);}
 535   
 536       tree->GetEntry(l);
 537       seltree->GetEntry(l);
 538       //  xsectree->GetEntry(l);
 539       //  MCLumi_Run = MCLumi_Run*factor;         
 540       weight = GetEventWeight(lumi,MCLumi_Run,P_e,P_p)*EventWeight;
 541       
 542       //cout << "EventWeight = "<<EventWeight<<"\t";
 543       //cout << "GetEventWeight = "<<(lumi/MCLumi_Run)*0.25<<"\t";
 544       //printf("w   %.3e\n",weight);
 545       //  Pe = 0.0;
 546       //  Pp = 0.0;
 547      
 548      
 549       int _bin = h_bg_dsigmadE->GetXaxis()->FindBin(SelPhoton_energy[0]); // dsigma/dE in fb
 550       xsec_bg = h_bg_dsigmadE->GetBinContent(_bin);
 551       //cout << "bg " << xsec_bg << endl;
 552       _bin = h_WIMPdsigmadE[0][0][0]->GetXaxis()->FindBin(SelPhoton_energy[0]);
 553 
 554   
 555       for (int m(0); m < 9 ; m++){
 556 
 557 	if (m == 0) {
 558   	  Pe = 0.0;
 559   	  Pp = 0.0;
 560   	}
 561   	if (m == 1) {
 562   	  Pe = 0.8;
 563   	  Pp = 0.0;
 564   	}
 565   	if (m == 2) {
 566   	  Pe = 0.8;
 567   	  Pp = 0.3;//Pp = -0.3
 568   	}
 569   	if (m == 3) {
 570   	  Pe = 0.8;
 571   	  Pp = 0.6;//Pp = -0.6
 572   	}
 573   	if (m == 4) {
 574   	  Pe = -0.8;
 575   	  Pp = 0.0;
 576   	}
 577   	if (m == 5) {
 578   	  Pe = -0.8;
 579   	  Pp = 0.3;
 580   	}
 581   	if (m == 6) {
 582   	  Pe = -0.8;
 583   	  Pp = 0.6;
 584   	}
 585         if (m == 7) {
 586           Pe = 0.8;
 587           Pp = -0.3;
 588        }
 589         if (m == 8) {
 590           Pe = 0.8;
 591           Pp = -0.6;
 592        }
 593     
 594   	for (int i(0); i < 250; i++){
 595   	  for (int b(0);  b< 1; b++){// b == l ( for _lambda)
 596   	         int s=0;
 597 		 xsec_sig = h_WIMPdsigmadE[i][s][b]->GetBinContent(_bin);
 598 	      //if(i==0)
 599 	      //cout << "sig  " << xsec_sig << endl;
 600   	     
 601 
 602   		if (k == 0){
 603 		//Vector operator
 604 		  _sigma_eRpL = 1.0;
 605   		  _sigma_eRpR = 0.0;
 606   		  _sigma_eLpL = 0.0;
 607   		  _sigma_eLpR = 1.0;
 608 		  		
 609 		}else
 610   		if (k == 1||k == 2){
 611 		//Axial-Vector operator And Scalar s-channel
 612 		  _sigma_eRpL = 0.0;
 613   		  _sigma_eRpR = 1.0;
 614   		  _sigma_eLpL = 1.0;
 615   		  _sigma_eLpR = 0.0;
 616   		}
 617   		
 618 		      
 619 		  		
 620 		_BR_eff = _BR * Sigma_PolFactor(_sigma_eRpL, _sigma_eRpR, _sigma_eLpL,_sigma_eLpR,Pe,Pp); 
 621 		
 622 // 		if((_BR_eff!=1) && (i==10) ){
 623 // 		cout<<_sigma_eRpL << _sigma_eRpR<< _sigma_eLpL<<_sigma_eLpR<<endl;
 624 // 		cout<<" _BR  = "<<_BR;
 625 // 		cout<<" _BR_eff = "<<_BR_eff;
 626 // 		cout<<"  Pe = "<<Pe;
 627 // 		cout<<"  Pp = "<<Pp<<endl;
 628 // 		cout<<"  xsec_bg  = "<<xsec_bg<<endl;}
 629 		
 630   		if (xsec_bg > 0.0){
 631 		  double weight2 = weight * _BR_eff * xsec_sig/xsec_bg;
 632 		  //double weight2 = weight * xsec_sig/xsec_bg;
 633 		  // cout << "w  " << weight2 << endl;
 634 		  h_weight[m]->Fill(TMath::Cos(Theta_prime),weight2);
 635 		  h_array[m][i][b][0]->Fill(SelPhoton_energy[0], weight2);
 636 		}
 637   	      
 638   	    }
 639 	 }	     		      
 640       }
 641     }
 642     delete f; 
 643   }
 644 
 645   cout << "asdad\n";
 646   // // ###################################################################
 647   // // ###############            END           ##########################
 648   // // ###################################################################
 649 
 650   // for (int l(1); l < 101; l++);
 651 
 652   // h_weight->DrawCopy();
 653  
 654 
 655   for (int m(0); m < 9 ; m++){ 
 656 
 657     f_sig[m]->cd();
 658     h_weight[m]->Write();
 659    
 660     for (int i(0); i < 250; i++){
 661       for (int l(0); l < 1; l++){ // Lambda
 662 	
 663 	  k=0;
 664 	    h_array[m][i][l][k]->Write();
 665 	  	
 666       }
 667     }
 668      cout << "asdad\n";
 669     f_sig[m]->Close();
 670   }
 671  
 672   return;
 673 }
 674   
 675 
 676 
 677 
 678 
 679 
 680 	       // _bin = h_WIMPdsigmadE[i][s][j]->GetXaxis()->FindBin(SelPhoton_energy[0]);
 681 	       // xsec_sig =_BR_eff * h_WIMPdsigmadE[i][s][j]->GetBinContent(_bin);
 682 
 683     // s_prime = 500.0*500.0;
 684     // x = SelPhoton_energy[0]*2.0/sqrt(s_prime);
 685     // Theta_prime = TMath::ACos(SelPhoton_cosTheta[0]);
 686     // differentialxsection_(&s_prime,&x, &Theta_prime,&Pe,&Pp,&xsec_bg);
 687 
 688   //----------------------------------------------------------------------------
 689 
 690 
 691 
 692 
 693 
 694 
 695 
 696 
 697 
 698 
 699 
 700 
 701 
 702 	
 703     
 704 // // ###################################################################
 705 // // ################ CREATE SIGNAL TEMPLATES ##########################
 706 // // ################### METHOD 2, USE 1D BG XSEC ######################
 707  
 708 
 709 
 710 // ifstream txtfile2("../steeringfiles/ILD00_Signal_RootFiles.txt");
 711 // 
 712 // while (!txtfile2.eof()){
 713 //   txtfile2.getline(name,256,'\n');
 714 //   
 715 //   printf("%s\n",name);
 716 //   f = new TFile(name,"OPEN");
 717 //   if (f->IsZombie()) continue;
 718 //   
 719 //   // ********** GET BACKGROUND DATA TREES
 720 //   
 721 //   TTree* tree = (TTree*)f->Get("Tree");
 722 //   TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 723 //   TTree* xsectree = (TTree*)f->Get("XSecTree");
 724 //   LoadTree(tree);
 725 //   LoadSelPhotonTree(seltree);
 726 //   LoadXSecTree(xsectree);
 727 //   
 728 //   // #########  LOOP OVER ENTRIES
 729 //   int nentries = seltree->GetEntries();
 730 //   // nentries = 1e3;
 731 //
 732 //   for (int k(0); k < nentries; k++){
 733 //     tree->GetEntry(k);
 734 //     seltree->GetEntry(k);
 735 //     xsectree->GetEntry(k);
 736 //     if (k%10000 == 0) printf("event %i\n",k);
 737 // 
 738 //     int bin = h_bg_dsigmadE->GetXaxis()->FindBin(SelPhoton_energy[0]);
 739 //     weight = h_bg_dsigmadE->GetBinContent(bin);
 740 //     if ( weight > 0.0){
 741 //	weight = 1/weight;
 742 //	weight *= h_WIMPdsigmadE->GetBinContent(bin);
 743 //	h_weight->Fill(TMath::Cos(Theta_prime),weight);
 744 //	h_array[20][2][1][0]->Fill(SelPhoton_energy[0], weight);
 745 //     }
 746 //   }
 747 //   delete f;
 748 // } 
 749 //
 750 
 751 // // ###################################################################
 752 // // ################ CREATE SIGNAL TEMPLATES ##########################
 753 // // ################### METHOD 3, USE 2D BG XSEC ######################
 754  
 755 // ifstream txtfile2("../steeringfiles/ILD00_Signal_RootFiles.txt");
 756 // 
 757 //
 758 //
 759 // while (!txtfile2.eof()){
 760 //   txtfile2.getline(name,256,'\n');
 761 //   
 762 //   printf("%s\n",name);
 763 //   f = new TFile(name,"OPEN");
 764 //   if (f->IsZombie()) continue;
 765 //   
 766 //   // ********** GET BACKGROUND DATA TREES
 767 //   
 768 //   TTree* tree = (TTree*)f->Get("Tree");
 769 //   TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 770 //   TTree* xsectree = (TTree*)f->Get("XSecTree");
 771 //   LoadTree(tree);
 772 //   LoadSelPhotonTree(seltree);
 773 //   LoadXSecTree(xsectree);
 774 //   
 775 //   // #########  LOOP OVER ENTRIES
 776 //   int nentries = seltree->GetEntries();
 777 //   // nentries = 1e3;
 778 //
 779 //   for (int k(0); k < nentries; k++){
 780 //     tree->GetEntry(k);
 781 //     seltree->GetEntry(k);
 782 //     xsectree->GetEntry(k);
 783 //     if (k%10000 == 0) printf("event %i\n",k);
 784 // 
 785 //     int bin_x = h_bg_2d->GetXaxis()->FindBin(E_prime);
 786 //     int bin_c = h_bg_2d->GetYaxis()->FindBin(TMath::Cos(Theta_prime));
 787 //
 788 //     weight = h_bg_2d->GetBinContent(bin_x,bin_c);
 789 //
 790 //     if ( weight > 0.0){
 791 //	weight = 1/weight;
 792 //	xsec_sig = DiffWIMPCrossSection(sqrt(s_prime), E_prime, Theta_prime);
 793 //	weight *= xsec_sig;
 794 //	h_weight->Fill(TMath::Cos(Theta_prime),weight);
 795 //	h_array[20][2][1][0]->Fill(SelPhoton_energy[0], weight);
 796 //     }
 797 //   }
 798 //   delete f;
 799 // } 
 800 
 801 
 802 // // ###################################################################
 803 // // ################ CREATE SIGNAL TEMPLATES ##########################
 804 // // ################### METHOD 4, USE STORED BG XSEC ##################
 805 
 806 // ifstream txtfile2("../steeringfiles/ILD00_Signal_RootFiles.txt");
 807 //
 808 // while (!txtfile2.eof()){
 809 //   txtfile2.getline(name,256,'\n');
 810 //   
 811 //   printf("%s\n",name);
 812 //   f = new TFile(name,"OPEN");
 813 //   if (f->IsZombie()) continue;
 814 //   
 815 //   // ********** GET BACKGROUND DATA TREES
 816 //   
 817 //   TTree* tree = (TTree*)f->Get("Tree");
 818 //   TTree* seltree = (TTree*)f->Get("SelPhotonTree");
 819 //   TTree* xsectree = (TTree*)f->Get("XSecTree");
 820 //   LoadTree(tree);
 821 //   LoadSelPhotonTree(seltree);
 822 //   LoadXSecTree(xsectree);
 823 //
 824 //   //########################################    
 825 //   // #########  LOOP OVER ENTRIES
 826 //   
 827 //   int nentries = seltree->GetEntries();
 828 //   //nentries = 1e3;
 829 //   
 830 //   for (int k(0); k < nentries; k++){
 831 //     tree->GetEntry(k);
 832 //     seltree->GetEntry(k);
 833 //     xsectree->GetEntry(k);
 834 //     if (k%10000 == 0) printf("event %i\n",k);
 835 //       
 836 //     //for (int i(0); i < 40; i++){
 837 //     //	for (int s(2); s < 3; s++){
 838 //     //	  for ( int j(1); j < 2; j++){
 839 //     //	    _mass = (double)i + 100.0;
 840 //     //	    _spin = (double)s/2.0;
 841 //     //	    _annihilator = (double)j;
 842 //     //	    
 843 //     //	    if (j == 1) _sigma_an = 6.0;
 844 //     //	    else   _sigma_an = 0.85; 
 845 //     //	   
 846 //     //	    _BR_eff = _BR * GetWIMP_PolFactor(_kappa_eRpL, _kappa_eRpR, _kappa_eLpL,_kappa_eLpR,Pe,Pp);
 847 //     //	    // _BR_eff = _BR * GetWIMP_PolFactor(_kappa_eRpL, _kappa_eRpR, _kappa_eLpL,_kappa_eLpR,P_e,P_p);
 848 //     //	    
 849 //     //	    double weight = GetEventWeight(lumi,MCLumi_Run,P_e,P_p);
 850 //     //	    weight *= GetSignalWeight(sqrt(s_prime),E_prime,Theta_prime,XSec);
 851 //     //	    h_weight->Fill(TMath::Cos(Theta_prime),weight);	    
 852 //     //	    h_array[i][s][j][0]->Fill(SelPhoton_energy[0], weight);
 853 //     //	  }
 854 //     //	}
 855 //     //}
 856 //
 857 //     _BR_eff = _BR * GetWIMP_PolFactor(_kappa_eRpL, _kappa_eRpR, _kappa_eLpL,_kappa_eLpR,P_e,P_p);
 858 //     
 859 //     double weight = GetEventWeight(lumi,MCLumi_Run,P_e,P_p);
 860 //     weight *= GetSignalWeight(sqrt(s_prime),E_prime,Theta_prime,XSec);
 861 //     h_weight->Fill(TMath::Cos(Theta_prime),weight);
 862 //     h_array[20][2][1][0]->Fill(SelPhoton_energy[0], weight);
 863 //   }
 864 //   delete f;
 865 // }   
 866  

Attached Files

To refer to attachments on a page, use attachment:filename, as shown below in the list of files. Do NOT use the URL of the [get] link, since this is subject to change and can break easily.
  • [get | view] (2012-01-30 13:28:08, 4.3 KB) [[attachment:CalibrateData.C]]
  • [get | view] (2012-01-30 13:27:40, 12.0 KB) [[attachment:CalibrationFunctions.C]]
  • [get | view] (2014-06-20 09:25:03, 5.2 KB) [[attachment:CreateDataBackground.C]]
  • [get | view] (2014-06-20 09:21:28, 22.3 KB) [[attachment:CreateDataForSignalAndBackground.C]]
  • [get | view] (2014-06-20 09:24:20, 0.2 KB) [[attachment:CreateDataForSignalAndBackground.py]]
  • [get | view] (2012-01-30 12:32:43, 12.0 KB) [[attachment:FindMergeAngle.C]]
  • [get | view] (2012-01-30 13:22:33, 9.8 KB) [[attachment:PreprocessData.C]]
  • [get | view] (2012-01-30 13:51:59, 16.8 KB) [[attachment:Selection.C]]
  • [get | view] (2014-06-20 09:26:52, 13.2 KB) [[attachment:TSysLimit_searches.C]]
  • [get | view] (2014-06-20 09:28:25, 14.9 KB) [[attachment:TSysLimit_searches_onebin.C]]
  • [get | view] (2014-06-20 09:25:07, 12.9 KB) [[attachment:WIMP_functions.C]]
  • [get | view] (2014-06-20 12:52:10, 1.3 KB) [[attachment:runLambdaSearch.sh]]
  • [get | view] (2014-06-20 12:52:16, 0.2 KB) [[attachment:submitter.sh]]
 All files | Selected Files: delete move to page copy to page

You are not allowed to attach a file to this page.