46         cout << 
" lmd hit noise studies " << endl;
 
   47         TCanvas* canvas = 
new TCanvas(
"canvas", 
"canvas", 900, 600);
 
   52         tMC.Add(
"./Lumi_digi_0.root");
 
   54         TChain thits(
"pndsim");
 
   55         thits.Add(
"./Lumi_recoMerged_0.root");
 
   57         TChain ttrk(
"pndsim");
 
   58         ttrk.Add(
"./Lumi_TrackNotFiltered_0.root");
 
   60         TChain ttrk_filt(
"pndsim");
 
   61         ttrk_filt.Add(
"./Lumi_Track_0.root");
 
   63         bool usefiltered(
true);
 
   64         TChain ttrk_prop(
"pndsim");
 
   66                 ttrk_prop.Add(
"./Lumi_GeaneFiltered_0.root");
 
   68                 ttrk_prop.Add(
"./Lumi_Geane_0.root");
 
   75         TClonesArray* digis = 
new TClonesArray(
"PndSdsDigiPixel");
 
   76         tMC.SetBranchAddress(
"LMDPixelDigis", &digis); 
 
   78         TClonesArray* 
hits = 
new TClonesArray(
"PndSdsMergedHit");
 
   79         thits.SetBranchAddress(
"LMDHitsMerged", &hits); 
 
   81         TClonesArray* tracks = 
new TClonesArray(
"PndTrack");
 
   82         ttrk.SetBranchAddress(
"LMDPndTrack", &tracks);
 
   84         TClonesArray* tracks_filt = 
new TClonesArray(
"PndTrack");
 
   85         ttrk_filt.SetBranchAddress(
"LMDPndTrack", &tracks_filt);
 
   87         TClonesArray* tracks_prop = 
new TClonesArray(
"FairTrackParH");
 
   89                 ttrk_prop.SetBranchAddress(
"LMDCleanTrack", &tracks_prop);
 
   90                 cout << 
" using geanie filtered tracks " << endl;
 
   92                 ttrk_prop.SetBranchAddress(
"GeaneTrackFinal", &tracks_prop);
 
  122         TH2D* hist_track_dir = 
new TH2D(
"hist_track_dir", 
"track dir", 500, -200, 200, 500, -200, 200);
 
  123         TH2D* hist_track_dir_filt = 
new TH2D(
"hist_track_dir_filt", 
"track dir filtered", 500, -200, 200, 500, -200, 200);
 
  124         TH2D* hist_track_dir_prop = 
new TH2D(
"hist_track_dir_prop", 
"track dir propagated", 50, -20, 20, 50, -20, 20);
 
  125         TH1D* hist_track_theta = 
new TH1D(
"hist_track_theta", 
"#theta distr", 50, 0, 15);
 
  127         TH2D* hist_hit_mult = 
new TH2D(
"hist_hit_mult", 
"hit multiplicity", 1000, 0, 1000, 4, 0, 4);
 
  128         TH1D* hist_trk_mult = 
new TH1D(
"hist_trk_mult", 
"track multiplicity", 100, 0, 100);
 
  129         TH2D* hist_hit_distr = 
new TH2D(
"hist_hit_distr", 
"hit distribution", 100, -10, 10, 100, -10, 10);
 
  131         int nEvents = tMC.GetEntries();
 
  132         int nEvents_counted(0);
 
  135         int ntrackstotal_filt(0);
 
  136         int ntrackstotal_prop(0);
 
  137         int ntrackstotal_rel(0);
 
  139         cout << 
" reading " << nEvents << 
" Events " << endl;
 
  143         for (Int_t j = 0; j < 
nEvents ; j++) {
 
  148                 ttrk_filt.GetEntry(j);
 
  149                 ttrk_prop.GetEntry(j);
 
  151                 ttrk_prop.GetEntry(j);
 
  153                 int nTracks = tracks->GetEntriesFast();
 
  154                 hist_trk_mult->Fill(nTracks);
 
  156                 for (
auto iTrack = 0; iTrack < nTracks; iTrack++){
 
  163                                 hist_track_dir->Fill(pxpz, pypz);
 
  164                                 double theta2 = pxpz*pxpz + pypz*pypz;
 
  165                                 if ((2*2 < theta2 && theta2 < 12*12)){
 
  172                 nTracks = tracks_filt->GetEntriesFast();
 
  174                 for (
auto iTrack = 0; iTrack < nTracks; iTrack++){
 
  181                                 hist_track_dir_filt->Fill(pxpz, pypz);
 
  185                 nTracks = tracks_prop->GetEntriesFast();
 
  187                 for (
auto iTrack = 0; iTrack < nTracks; iTrack++){
 
  188                         FairTrackParH* track = (FairTrackParH*) tracks_prop->At(iTrack);
 
  190                           if(track->GetLambda()!=0){ 
 
  192                             double pxpz = track->GetPx()/track->GetPz()*1e3;
 
  193                             double pypz = track->GetPy()/track->GetPz()*1e3;
 
  194                             hist_track_dir_prop->Fill(pxpz, pypz);
 
  195                             double theta = 
sqrt(pxpz*pxpz+pypz*pypz);
 
  196                             hist_track_theta->Fill(theta);
 
  201                 int nHits = digis->GetEntriesFast();
 
  204                 int ihits[4] = {0,0,0,0};
 
  205                 for (Int_t iHit = 0; iHit < 
nHits; iHit++) {
 
  214                                 int ihalf, iplane, imodule, iside, idie, isens;
 
  216                                 lmddim.
Get_sensor_by_id(sensid, ihalf, iplane, imodule, iside, idie, isens);
 
  253                 hist_hit_mult->Fill(ihits[0], 0);
 
  254                 hist_hit_mult->Fill(ihits[1], 1);
 
  255                 hist_hit_mult->Fill(ihits[2], 2);
 
  256                 hist_hit_mult->Fill(ihits[3], 3);
 
  258                 nHits = hits->GetEntriesFast();
 
  259                 for (
auto iHit = 0; iHit < 
nHits; iHit++){
 
  262                         int ihalf, iplane, imodule, iside, idie, isensor;
 
  263                         lmddim.
Get_sensor_by_id(sensID, ihalf, iplane, imodule, iside, idie, isensor);
 
  266                         hist_hit_distr->Fill(hit_pos_lmd.X(), hit_pos_lmd.Y());
 
  283         cout << nEvents_counted << 
" events " << endl;
 
  284         cout << 
" in " << nhitstotal << 
" hits found " << real_hits << 
" real hits and " << noise_hits << 
" noise hits " << endl;
 
  285         cout << 
" the ratio over nevents is " << (double) noise_hits / (
double) nEvents_counted << endl;
 
  286         cout << " found " << ntrackstotal << " tracks " << (
double) ntrackstotal/(
double) nEvents_counted << endl;
 
  287         cout << " found " << ntrackstotal_rel << " tracks relevant 
for analysis " << (
double) ntrackstotal_rel/(
double) nEvents_counted << endl;
 
  288         cout << " found " << ntrackstotal_filt << " filtered tracks " << (
double) ntrackstotal_filt/(
double) nEvents_counted << endl;
 
  289         cout << " found " << ntrackstotal_prop << " propagated tracks " << (
double) ntrackstotal_prop/(
double) nEvents_counted << endl;
 
  292         hist_track_dir->
Draw("COLZ");
 
  294         hist_track_dir->SetYTitle("p_{
y}/p_{
z}x10^{-3}
"); 
  297         hist_hit_distr->Draw("COLZ
"); 
  298         hist_hit_distr->SetXTitle("X [
cm]
"); 
  299         hist_hit_distr->SetYTitle("Y [
cm]
"); 
  301         //hist_track_dir_filt->Draw("COLZ
"); 
  302         //hist_track_dir_filt->SetXTitle("p_{
x}/p_{
z}x10^{-3}
"); 
  303         //hist_track_dir_filt->SetYTitle("p_{
y}/p_{
z}x10^{-3}
"); 
  306         hist_trk_mult->Draw(); 
  307         hist_trk_mult->SetXTitle("# tracks / 
event"); 
  310         hist_hit_mult->Draw("colz
"); 
  311         hist_hit_mult->SetXTitle("# hits / 
event"); 
  312         hist_hit_mult->SetYTitle("plane
"); 
  315         hist_track_dir_prop->Draw("COLZ
"); 
  316         hist_track_dir_prop->SetXTitle("p_{
x}/p_{
z}x10^{-3}
"); 
  317         hist_track_dir_prop->SetYTitle("p_{
y}/p_{
z}x10^{-3}
"); 
  320         hist_track_theta->Draw(); 
  321         hist_track_theta->SetXTitle("#Theta [mrad]
"); 
  323         plane3_rate->Scale(1./tsim*1e-3); // in kHz 
  324         //plane3_rate->SetMaximum(150); 
  325         cout << " total rate of 
" << nhitstotal << " hits 
is " << nhitstotal/tsim*1e-6 << "MHz
" << endl; 
  326         plane3_rad->Scale(1./tsim*60.*60.*24.*365.*0.5/(2.33e-3*(2*2*50.e-4))); // 50% duty cycle over a year and normalized to the mass in kg of one sensor 
  328         plane3_Edep->Draw("colz
"); 
  329         //hist_eloss_total->Draw(); 
  331         plane3_hits->Draw("colz
"); 
  332         canvas->Print("MC_hit_distr_1_5_plane_3_side_0.pdf
"); 
  333         canvas->Print("MC_hit_distr_1_5_plane_3_side_0.root
"); 
  335         canvas->Print("MC_IEL_distrib_1_5.pdf
"); 
  336         canvas->Print("MC_IEL_distrib_1_5.root
"); 
  338         plane3_rate->Draw("colz
"); 
  339         cout << " max rate 
is " << plane3_rate->GetBinContent(plane3_rate->GetMaximumBin()) << " kHz
" << endl; 
  340         canvas->Print("MC_hit_rate_1_5_plane_3_side_0.pdf
"); 
  341         canvas->Print("MC_hit_rate_1_5_plane_3_side_0.root
"); 
  343         plane3_rad->Draw("colz
"); 
  344         cout << " max dose 
is " << plane3_rad->GetBinContent(plane3_rad->GetMaximumBin()) << " Gy
" << endl; 
  345         canvas->Print("MC_IEL_dose_1_5_plane_3_side_0.pdf
"); 
  346         canvas->Print("MC_IEL_does_1_5_plane_3_side_0.root
"); 
Int_t GetSensorID() const 
TVector3 GetPosition() const 
friend F32vec4 sqrt(const F32vec4 &a)
TVector3 Transform_global_to_lmd_local(const TVector3 &point, bool isvector=false, bool aligned=true)
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
hist4 SetXTitle("p/(GeV/c)")
for(int j=0;j< ncounts;j++)
tree Draw("fELoss:TMath::Sqrt(fPx_out*fPx_out+ fPy_out*fPy_out+ fPz_out*fPz_out)", ppos &&"fELoss < 0.04")
void Get_sensor_by_id(const int sensor_id, int &ihalf, int &iplane, int &imodule, int &iside, int &idie, int &isensor)
int hit(Int_t nEvents=0, TString inFile="sim.root", TString parFile="par.root", TString inDigi="digi.root", TString outFile="hit.root", Int_t timeBased=0)
static PndLmdDim & Get_instance()
Data class to store the digi output of a pixel module. 
Int_t GetSensorID() const 
void DrawProgressBar(int len, double percent)
FairTrackParP GetParamFirst()