ROOT logo
/**
 * @defgroup pwglf_forward_scripts_corr Correction scripts 
 * 
 * Scripts and functions to handle corrections 
 *
 * @ingroup pwglf_forward_scripts
 */

//____________________________________________________________________
/** 
 * Compare secondary maps  
 * 
 * @param fn1  File 1
 * @param fn2  File 2
 * @param n1   Name 1
 * @param n2   Name 2
 * @param load 
 *
 * @ingroup pwglf_forward_scripts_corr
 */
void
CompareCentralSecMaps(const char* fn1,   const char* fn2, 
		      const char* n1=0,  const char* n2=0,
		      bool load=true)
{

  // --- Load Utilities ----------------------------------------------
  if (load) {
    gROOT->Macro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/LoadLibs.C");
    //   gROOT->LoadMacro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/CompareCorrs.C");
  }

  // --- Get Objects -------------------------------------------------
  TObject* o1 = 0;
  TObject* o2 = 0;
  TFile file1(fn1);
  TFile file2(fn2);
  
  AliCentralMultiplicityTask task("tmp");
  
  o1 = file1.Get(task.GetManager().GetSecMapName());
  o2 = file2.Get(task.GetManager().GetSecMapName());
  if (!o1 || !o2) return; 
  AliCentralCorrSecondaryMap* obj1 = static_cast<AliCentralCorrSecondaryMap*>(o1);
  AliCentralCorrSecondaryMap* obj2 = static_cast<AliCentralCorrSecondaryMap*>(o2);
  UShort_t nVtx = obj1->GetVertexAxis().GetNbins();

  // --- Make canvas -------------------------------------------------
  // Canvas* c = new Canvas("CentralSecMapComparison", "Ratio of central secondary maps", n1, n2);
  // c->Open();

  // --- Loop over the data ------------------------------------------
  //for (UShort_t d = 1; d <= 3; d++) { 
  //  UShort_t nR = (d == 1 ? 1 : 2);
  //  for (UShort_t q = 0; q < nR; q++) { 
  //   Char_t   r  = (q == 0 ? 'I' : 'O');
  //    UShort_t nS = (q == 0 ?  20 :  40);

      // --- Make 2D ratios ------------------------------------------
      //c->Clear(nVtx, d, r);
  TCanvas* c2D = new TCanvas("2Dcomparison","2Dcomparison",800,1000);
  TCanvas* c1D = new TCanvas("1Dcomparison","1Dcomparison",800,1000);
  c2D->Divide(2,5);
  c1D->Divide(2,5);
  c2D->cd();
  
  TList hists;
  for (UShort_t v=1; v <= nVtx; v++) { 
    //	TVirtualPad* p = c->cd(v);
    c2D->cd(v);
    TH2* h1 = obj1->GetCorrection( v);
    TH2* h2 = obj2->GetCorrection( v);
    //std::cout<<h1<<"  "<<h2<<std::endl;
    /*	if (!h1) { 
	Error("CompareSecMaps", 
	"Map for SPD, vtxbin %3d not found in first", 
	v);
	continue;
	}
	if (!h2) { 
	Error("CompareSecMaps", 
	"Map for SPD, vtxbin %3d not found in second", 
	v);
	continue;
	}
    */
    Double_t vl    = obj1->GetVertexAxis().GetBinLowEdge(v);
    Double_t vh    = obj1->GetVertexAxis().GetBinUpEdge(v);
    TH2*     ratio = 
      static_cast<TH2*>(h1->Clone(Form("tmpSPD_%3d",v)));
    ratio->SetName(Form("SPD_vtx%03d_ratio",  v));
    ratio->SetTitle(Form("%+5.1f<v_{z}<%-+5.1f", vl, vh));
    ratio->Divide(h2);
    ratio->SetStats(0);
    ratio->SetDirectory(0);
    ratio->SetZTitle("ratio");
    
	//if (ratio->GetMaximum()-ratio->GetMinimum() > 10) 
	//  p->SetLogz();

    ratio->DrawCopy("colz");
    hists.AddAt(ratio, v-1);
    
  }
      // c->Print(d, r);
      
      // --- Make 1D profiles ----------------------------------------
      //c->Clear(nVtx, d, r);
    
  c1D->cd();
  for (UShort_t v=1; v <= nVtx; v++) { 
    //c->cd(v);
    c1D->cd(v);
    TH2* hist = static_cast<TH2*>(hists.At(v-1));
    TH1D* prof = hist->ProjectionX();
    prof->Clear();
    for(Int_t i=1; i<=hist->GetNbinsX();i++) {
      
      Float_t sum = 0;
      Float_t error = 0;
      Int_t   n   = 0;
      for(Int_t j=1; j<=hist->GetNbinsY();j++) {
	if(hist->GetBinContent(i,j) > 0) {
	  sum = sum+hist->GetBinContent(i,j);
	  error = error + TMath::Power(hist->GetBinError(i,j),2);
	  n++;
	}
	
      }
      if(n>0) {
	sum = sum/(Float_t)n;
	error = TMath::Sqrt(error) / (Float_t)n;
	prof->SetBinContent(i,sum);
	prof->SetBinError(i,error);
      }
    }
    
    //    prof->Scale(0.05);
    prof->SetStats(0);
    prof->SetMinimum(0.8);
    prof->SetMaximum(1.2);
    
    
    prof->Fit("pol0","Q");
    prof->DrawCopy();
    
    TF1* f = prof->GetFunction("pol0");
    
    TLatex* l = new TLatex(0.5, 0.4, Form("A = %f #pm %f", 
					  f->GetParameter(0), 
					  f->GetParError(0)));
    l->SetTextAlign(22);
    l->SetNDC();
    l->Draw();
    l->DrawLatex(0.5, 0.3, Form("#chi^2/NDF = %f / %d = %f", 
    			f->GetChisquare(), 
    			f->GetNDF(), 
    			f->GetChisquare() / f->GetNDF()));
    Double_t dist = TMath::Abs(1 - f->GetParameter(0));
    l->DrawLatex(0.5, 0.35, Form("|1 - A| = %f %s #deltaA", 
	 dist, dist <= f->GetParError(0) ? 
		 "#leq" : ">")); 
    
    TLine* l1 = new TLine(-4, 1, 6, 1);
    l1->SetLineColor(kRed);
    l1->SetLineStyle(2);
    l1->Draw(); 
    
  }
  
  //c->Print(d, r, "profiles");
  c2D->Print("comparisonSPD.pdf(");
  c1D->Print("comparisonSPD.pdf)");
  

  // --- Close stuff -------------------------------------------------
  //c->Close();
  // file1->Close();
  // file2->Close();
}

  
//____________________________________________________________________
//
// EOF
// 
 CompareCentralSecMaps.C:1
 CompareCentralSecMaps.C:2
 CompareCentralSecMaps.C:3
 CompareCentralSecMaps.C:4
 CompareCentralSecMaps.C:5
 CompareCentralSecMaps.C:6
 CompareCentralSecMaps.C:7
 CompareCentralSecMaps.C:8
 CompareCentralSecMaps.C:9
 CompareCentralSecMaps.C:10
 CompareCentralSecMaps.C:11
 CompareCentralSecMaps.C:12
 CompareCentralSecMaps.C:13
 CompareCentralSecMaps.C:14
 CompareCentralSecMaps.C:15
 CompareCentralSecMaps.C:16
 CompareCentralSecMaps.C:17
 CompareCentralSecMaps.C:18
 CompareCentralSecMaps.C:19
 CompareCentralSecMaps.C:20
 CompareCentralSecMaps.C:21
 CompareCentralSecMaps.C:22
 CompareCentralSecMaps.C:23
 CompareCentralSecMaps.C:24
 CompareCentralSecMaps.C:25
 CompareCentralSecMaps.C:26
 CompareCentralSecMaps.C:27
 CompareCentralSecMaps.C:28
 CompareCentralSecMaps.C:29
 CompareCentralSecMaps.C:30
 CompareCentralSecMaps.C:31
 CompareCentralSecMaps.C:32
 CompareCentralSecMaps.C:33
 CompareCentralSecMaps.C:34
 CompareCentralSecMaps.C:35
 CompareCentralSecMaps.C:36
 CompareCentralSecMaps.C:37
 CompareCentralSecMaps.C:38
 CompareCentralSecMaps.C:39
 CompareCentralSecMaps.C:40
 CompareCentralSecMaps.C:41
 CompareCentralSecMaps.C:42
 CompareCentralSecMaps.C:43
 CompareCentralSecMaps.C:44
 CompareCentralSecMaps.C:45
 CompareCentralSecMaps.C:46
 CompareCentralSecMaps.C:47
 CompareCentralSecMaps.C:48
 CompareCentralSecMaps.C:49
 CompareCentralSecMaps.C:50
 CompareCentralSecMaps.C:51
 CompareCentralSecMaps.C:52
 CompareCentralSecMaps.C:53
 CompareCentralSecMaps.C:54
 CompareCentralSecMaps.C:55
 CompareCentralSecMaps.C:56
 CompareCentralSecMaps.C:57
 CompareCentralSecMaps.C:58
 CompareCentralSecMaps.C:59
 CompareCentralSecMaps.C:60
 CompareCentralSecMaps.C:61
 CompareCentralSecMaps.C:62
 CompareCentralSecMaps.C:63
 CompareCentralSecMaps.C:64
 CompareCentralSecMaps.C:65
 CompareCentralSecMaps.C:66
 CompareCentralSecMaps.C:67
 CompareCentralSecMaps.C:68
 CompareCentralSecMaps.C:69
 CompareCentralSecMaps.C:70
 CompareCentralSecMaps.C:71
 CompareCentralSecMaps.C:72
 CompareCentralSecMaps.C:73
 CompareCentralSecMaps.C:74
 CompareCentralSecMaps.C:75
 CompareCentralSecMaps.C:76
 CompareCentralSecMaps.C:77
 CompareCentralSecMaps.C:78
 CompareCentralSecMaps.C:79
 CompareCentralSecMaps.C:80
 CompareCentralSecMaps.C:81
 CompareCentralSecMaps.C:82
 CompareCentralSecMaps.C:83
 CompareCentralSecMaps.C:84
 CompareCentralSecMaps.C:85
 CompareCentralSecMaps.C:86
 CompareCentralSecMaps.C:87
 CompareCentralSecMaps.C:88
 CompareCentralSecMaps.C:89
 CompareCentralSecMaps.C:90
 CompareCentralSecMaps.C:91
 CompareCentralSecMaps.C:92
 CompareCentralSecMaps.C:93
 CompareCentralSecMaps.C:94
 CompareCentralSecMaps.C:95
 CompareCentralSecMaps.C:96
 CompareCentralSecMaps.C:97
 CompareCentralSecMaps.C:98
 CompareCentralSecMaps.C:99
 CompareCentralSecMaps.C:100
 CompareCentralSecMaps.C:101
 CompareCentralSecMaps.C:102
 CompareCentralSecMaps.C:103
 CompareCentralSecMaps.C:104
 CompareCentralSecMaps.C:105
 CompareCentralSecMaps.C:106
 CompareCentralSecMaps.C:107
 CompareCentralSecMaps.C:108
 CompareCentralSecMaps.C:109
 CompareCentralSecMaps.C:110
 CompareCentralSecMaps.C:111
 CompareCentralSecMaps.C:112
 CompareCentralSecMaps.C:113
 CompareCentralSecMaps.C:114
 CompareCentralSecMaps.C:115
 CompareCentralSecMaps.C:116
 CompareCentralSecMaps.C:117
 CompareCentralSecMaps.C:118
 CompareCentralSecMaps.C:119
 CompareCentralSecMaps.C:120
 CompareCentralSecMaps.C:121
 CompareCentralSecMaps.C:122
 CompareCentralSecMaps.C:123
 CompareCentralSecMaps.C:124
 CompareCentralSecMaps.C:125
 CompareCentralSecMaps.C:126
 CompareCentralSecMaps.C:127
 CompareCentralSecMaps.C:128
 CompareCentralSecMaps.C:129
 CompareCentralSecMaps.C:130
 CompareCentralSecMaps.C:131
 CompareCentralSecMaps.C:132
 CompareCentralSecMaps.C:133
 CompareCentralSecMaps.C:134
 CompareCentralSecMaps.C:135
 CompareCentralSecMaps.C:136
 CompareCentralSecMaps.C:137
 CompareCentralSecMaps.C:138
 CompareCentralSecMaps.C:139
 CompareCentralSecMaps.C:140
 CompareCentralSecMaps.C:141
 CompareCentralSecMaps.C:142
 CompareCentralSecMaps.C:143
 CompareCentralSecMaps.C:144
 CompareCentralSecMaps.C:145
 CompareCentralSecMaps.C:146
 CompareCentralSecMaps.C:147
 CompareCentralSecMaps.C:148
 CompareCentralSecMaps.C:149
 CompareCentralSecMaps.C:150
 CompareCentralSecMaps.C:151
 CompareCentralSecMaps.C:152
 CompareCentralSecMaps.C:153
 CompareCentralSecMaps.C:154
 CompareCentralSecMaps.C:155
 CompareCentralSecMaps.C:156
 CompareCentralSecMaps.C:157
 CompareCentralSecMaps.C:158
 CompareCentralSecMaps.C:159
 CompareCentralSecMaps.C:160
 CompareCentralSecMaps.C:161
 CompareCentralSecMaps.C:162
 CompareCentralSecMaps.C:163
 CompareCentralSecMaps.C:164
 CompareCentralSecMaps.C:165
 CompareCentralSecMaps.C:166
 CompareCentralSecMaps.C:167
 CompareCentralSecMaps.C:168
 CompareCentralSecMaps.C:169
 CompareCentralSecMaps.C:170
 CompareCentralSecMaps.C:171
 CompareCentralSecMaps.C:172
 CompareCentralSecMaps.C:173
 CompareCentralSecMaps.C:174
 CompareCentralSecMaps.C:175
 CompareCentralSecMaps.C:176
 CompareCentralSecMaps.C:177
 CompareCentralSecMaps.C:178
 CompareCentralSecMaps.C:179
 CompareCentralSecMaps.C:180
 CompareCentralSecMaps.C:181
 CompareCentralSecMaps.C:182
 CompareCentralSecMaps.C:183
 CompareCentralSecMaps.C:184
 CompareCentralSecMaps.C:185
 CompareCentralSecMaps.C:186
 CompareCentralSecMaps.C:187