musrsim/run/macros/td_mcp2_tof.C
2025-03-27 22:05:05 +01:00

106 lines
3.7 KiB
C

// To run root without logon screen and script use: root -l -n
// To run a named macro in root use:
//.L td_mcp2_tof.C
// td_mcp2_tof("data/sr1_1100.root")
#include <string>
#define NSIZE 100000
void td_mcp2_tof(const char* fname)
{
Double_t *ptofTD, *ptofMCP2, *pEventIDTD, *pEventIDMCP2;
Double_t tofTD[NSIZE], tofMCP2[NSIZE], EventIDTD[NSIZE], EventIDMCP2[NSIZE];
Long64_t i, j, k, nEventsTD, nEventsMCP2;
TTree *t1;
TH1D *htof, *tof_mcp2_plane, *tof_td;
TFile* f1 = new TFile(fname);
t1 = (TTree*) gDirectory->FindObjectAny("t1");
// muon time-of-flight to TD
t1->Draw("eventID:1000*save_time","save_detID==890&&((save_particleID==-13)||(save_particleID==-1313))"); //0, 2
// t1->Draw("eventID:1000*save_time","save_detID==890&&save_particleID==-1313"); //0, 2
nEventsTD = t1->GetSelectedRows();
ptofTD = t1->GetV2();
pEventIDTD = t1->GetV1();
for (i = 0; i<nEventsTD; i++){
EventIDTD[i] = pEventIDTD[i];
tofTD[i] = ptofTD[i];
}
// muon time-of-flight to MCP2, mu+ or Mu
t1->Draw("eventID:1000*save_time","save_detID==902&&((save_particleID==-13)||(save_particleID==-1313))&&sqrt(save_x**2+save_y**2)<=20"); //0, 2
// t1->Draw("eventID:1000*save_time","save_detID==902&&save_particleID==-1313"); //0, 2
nEventsMCP2 = t1->GetSelectedRows();
ptofMCP2 = t1->GetV2();
pEventIDMCP2 = t1->GetV1();
// cout << "nEventsMCP2 = " << nEventsMCP2 << endl;
for (i = 0; i<nEventsMCP2; i++){
EventIDMCP2[i] = pEventIDMCP2[i];
tofMCP2[i] = ptofMCP2[i];
}
// TOF TD-MCP2
htof = new TH1D("tof_TDMCP2", "Time-of-flight TD-MCP2", 1000, 0.25, 500.25);
k = 0;
for (i = 0; i<nEventsMCP2; i++){
for (j = k; j<nEventsTD; j++){
// cout << i << ", " << j << ", " << EventIDMCP2[i] << ", " << EventIDTD[j] << ", " << tofMCP2[i] << ", " << tofTD[j] << endl;
if (EventIDMCP2[i] == EventIDTD[j]){
htof->Fill(tofMCP2[i] - tofTD[j]);
k = j + 1;
break;
}
}
}
gStyle->SetPalette(1,0);
// gStyle->SetName(fname);
gStyle->SetOptTitle(1);
TCanvas* c1 = new TCanvas("Time of Flight","Time of Flight",60,40,600,800);
c1->Divide(1,3);
c1->SetGridx();
c1->SetGridy();
// muon time-of-flight to TD
c1->cd(1);
t1->Draw("1000*save_time >> tof_td(2000, 0.25, 1000.25)","save_detID==890&&save_particleID==-13"); //0, 2
tof_td = (TH1D*) gDirectory->Get("tof_td");
tof_td->SetTitle("Time-of-flight to TD");
tof_td->GetXaxis()->SetTitle("TD time of flight (ns)");
// tof_td->GetXaxis()->SetNdivisions(405);
tof_td->GetYaxis()->SetNdivisions(406);
// tof_td->GetXaxis()->SetTickLength(0.018);
tof_td->GetYaxis()->SetTickLength(0.018);
tof_td->Draw("F");
tof_td->SetFillStyle(1001);
tof_td->SetFillColor(kGreen-5);
tof_td->Draw();
//c1->SaveAs("data/sr_bspot_0.eps");
// muon time-of-flight to MCP2
c1->cd(2);
t1->Draw("1000*save_time >> tof_mcp2_plane(2000, 0.25, 1000.25)","save_detID==902&&((save_particleID==-13)||(save_particleID==-1313))"); //0, 2
tof_mcp2_plane = (TH1D*) gDirectory->Get("tof_mcp2_plane");
tof_mcp2_plane->SetTitle("Time-of-flight to plane of MCP2 (including muons missing MCP2)");
tof_mcp2_plane->GetXaxis()->SetTitle("MCP2 time of flight (ns)");
// tof_mcp2->GetXaxis()->SetNdivisions(405);
tof_mcp2_plane->GetYaxis()->SetNdivisions(406);
// tof_mcp2->GetXaxis()->SetTickLength(0.018);
tof_mcp2_plane->GetYaxis()->SetTickLength(0.018);
tof_mcp2_plane->Draw("F");
tof_mcp2_plane->SetFillStyle(1001);
tof_mcp2_plane->Draw();
c1->cd(3);
htof->SetTitle("Time-of-flight TD-MCP2");
htof->GetXaxis()->SetTitle("TD-MCP2 time of flight (ns)");
htof->GetYaxis()->SetNdivisions(406);
htof->GetYaxis()->SetTickLength(0.018);
htof->Draw("F");
htof->Draw();
}