gamma.cxx 3.75 KB
Newer Older
1 2 3 4 5
TChain* MakeChain1();
TChain* MakeChain2();
TChain* MakeChain();
TH2F*   MakeTH2();
TH2F*   GetTH2();
6
TGraph* graph = new TGraph(200);
7 8 9 10 11 12 13 14
double off;
double c_light=299.792458;//mm/ns
auto chain = MakeChain();
void process1bar(int b);
auto h = new TH2F("h","h",200,0,200,500,0,1000);
auto r = new TH2F("r","r",200,0,200,10000,17000,20000);
ofstream output("Calibration/Nebula/offset_gamma.txt");
////////////////////////////////////////////////////////////////////////////////
Adrien Matta's avatar
Adrien Matta committed
15
void gamma(){
16 17 18 19 20 21 22 23 24 25 26 27
  double beta=0.5504;
  // Distance from SBT(-7377.56) to target(3774.7 for FP + a value for each target)
  double D1=7377.56-3774.7+2.795;
  double D2=7377.56-3774.7+6.972;
  // time offset for each case
  double off1 = D1/(beta*c_light);
  double off2 = D2/(beta*c_light);
  off = (off1+off2)*0.5;
  cout << "Offset1 : " << off1 << endl;
  cout << "Offset2 : " << off2 << endl;
  cout << "Mean : " << off << endl;
  auto c = new TCanvas("tof","tof");
Adrien Matta's avatar
Adrien Matta committed
28

29 30 31
  chain->SetAlias("R","sqrt(Nebula.PosX*Nebula.PosX+Nebula.PosY*Nebula.PosY+(Nebula.PosZ+3774.7)*(Nebula.PosZ+3774.7))");
  h=GetTH2();
  h->Draw("colz");
Adrien Matta's avatar
Adrien Matta committed
32
  new TCanvas();
33 34 35 36
  for(unsigned int i = 0 ; i < 160 ; i++){
    process1bar(i); 
  }
  output.close();
37 38
  new TCanvas();
  graph->Draw("ap");
39 40 41 42
}

////////////////////////////////////////////////////////////////////////////////
void process1bar(int b){
Adrien Matta's avatar
Adrien Matta committed
43
  //new TCanvas();
44 45 46 47 48
  // Get the Radius for the distance to this barre
  auto r1 = r->ProjectionY(Form("h%d",b),b,b+1);  
  double R =  r1->GetBinCenter(r1->GetMaximumBin());

  auto h1 = h->ProjectionY(Form("h%d",b),b,b+1);
49
  h1->Rebin(8);
50
  double max = h1->GetBinCenter(h1->GetMaximumBin());
Adrien Matta's avatar
Adrien Matta committed
51
  //h1->Draw();
52
  auto f = new TF1("f","crystalball(0)+pol1(5)",max-50,max+50);
53 54
  f->SetParameter(0,h1->GetMaximum());
  f->SetParameter(1,max);
55 56 57 58
  f->SetParameter(2,35);
  f->SetParameter(3,0.1);
  f->SetParameter(4,1);

59 60
  h1->Fit(f,"R");
  
61 62 63 64 65 66
    // Vbad = R/(TOF) -> TOF/R = 1/Vbad
    // c= R/(TOF+X)   -> (TOF+X)/R = 1/c -> TOF/R+X/R=1/c
    // 1/Vbad +X/R = 1/c
    // X=R*(1/c-1/Vbad)
    
  double offset=R*(1/c_light-1/f->GetParameter(1)) ;
67 68
 // double offset=R*(1/c_light-1/max) ;
  
69
  cout <<f->GetParameter(1) << " " <<  offset << " " << R/(f->GetParameter(1)-offset) << endl;
70
  if(offset>0){
71
    output << "NEBULA_T_ID"  << b << " " << offset << endl; 
72 73
    graph->SetPoint(b,b,offset);
    }
74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99
}
////////////////////////////////////////////////////////////////////////////////
TH2F* GetTH2(){
  auto File= new TFile("Calibration/Nebula/hist_v.root");
  h = (TH2F*) File->FindObjectAny("h");

  File= new TFile("Calibration/Nebula/hist_r.root");
  r = (TH2F*) File->FindObjectAny("r");
  return h;
  }

////////////////////////////////////////////////////////////////////////////////
TH2F* MakeTH2(){
  TString cond=Form("(Nebula.TOF-%f)>20&&(Nebula.TOF-%f)<38",off,off);
  TString draw=Form("R/(Nebula.TOF-%f):Nebula.DetectorNumber>>h",off); 
  chain->Draw(draw,cond,"colz");
  h->SaveAs("Calibration/Nebula/hist_v.root");

  chain->Draw("R:Nebula.DetectorNumber>>r","","colz");
  r->SaveAs("Calibration/Nebula/hist_r.root");

  return h;
  }

////////////////////////////////////////////////////////////////////////////////
TChain* MakeChain(){
Adrien Matta's avatar
Adrien Matta committed
100
  auto chain = new TChain("PhysicsTree");
101 102 103
  chain->Add("root/analysis/gamma/run*.root");
  return chain;
}
Adrien Matta's avatar
Adrien Matta committed
104

105 106 107 108 109 110 111 112
////////////////////////////////////////////////////////////////////////////////
TChain* MakeChain1(){
  auto chain = new TChain("PhysicsTree");
  for(unsigned int i = 418 ; i < 442  ; i++){
    chain->Add(Form("root/analysis/gamma/run%d.root",i));
  }
  return chain;
}
Adrien Matta's avatar
Adrien Matta committed
113

114 115 116 117 118
////////////////////////////////////////////////////////////////////////////////
TChain* MakeChain2(){
  auto chain = new TChain("PhysicsTree");
  for(unsigned int i = 488 ; i < 498  ; i++){
    chain->Add(Form("root/analysis/gamma/run%d.root",i));
Adrien Matta's avatar
Adrien Matta committed
119
  }
120 121
  return chain;
}