/data3/calcul/jacquem/working_dir/Micromegas/micromegasFrameWork/src/analyse/TB2010/Testbeam_june_10/Complete_detector/Beam_profile_m2.C

Go to the documentation of this file.

{
  gROOT->Reset();
  gStyle->SetPalette(1);
  gSystem->Load("/LC/Detecteurs/MicroMegas/Offline/micromegasFrameWork/trunk/lib/libMicro.so");


  bool plot = 1;
  bool print = 1;
  bool signal_to_noise = 0;

  int nhit_cut = 1e9;

  TString file = "b0702_10_acq_HR2_02072010_1302_merged.root";
  //cout<<"HR2 file? ";
  //cin>>file;

  TString filename = "/lapp_data/LC/Detecteurs/MicroMegas/data/TB2010/SPS_H4_june_2010/Root_files/pro/"+file;

  int tested_chip = 124;


  float nTrees= 0;



  Long64_t t,t1,t2,t3,dt;
  int nchannel = 0;
  int hardid = 0;
  int softid = 0;
  int chipid = 0;
  int nhit = 0;
  int xpos = 0;  
  int ypos = 0;  
  int zpos = 0;  
  int adc = 0;  
  int nevent = 0;
  int nhit = 0;
  int nhit_intime = 0;
  int nhit_intime_tot = 0;

  TH1I * hdt = new TH1I("hdt","Time to trigger; #Deltat (clock cycle)",400,-10,390);    

  TH2I * hxy_m2 = new TH2I("hxy_m2","Raw profile;x (cm);y (cm)",96,0,96,96,0,96);
  TH2I * hxy_m2_cut = new TH2I("hxy_m2_cut","Profile after time cut;x (cm);y (cm)",96,0,96,96,0,96);

  TH1I * hxy_m2_cut_px = new TH1I("hxy_m2_cut_px","Profile after time cut;x (cm)",96,-0,96);
  TH1I * hxy_m2_cut_py = new TH1I("hxy_m2_cut_py","Profile after time cut;y (cm)",96,-0,96);

  TH1I * hnhit = new TH1I("hnhit","Number of hits",3200,0,3200);  
  TH1I * hnhit_intime = new TH1I("hnhit_intime","Number of his after time cut",3200,0,3200);  



  TFile *f = new TFile(filename);
  cout<<"reading File : "<<filename<<endl;

  TIter nextkey(f.GetListOfKeys());
  TKey *key;

  while (key = (TKey*)nextkey()) {
    
    TTree *tree = (TTree*)key->ReadObj();                
    MTRun* run = (MTRun*)tree->GetUserInfo()->FindObject("MTRun");
    //run->Info();
  nTrees++;
  float nEvent = tree.GetEntries();
  cout <<"nEvent="<<nEvent<<endl;

  MTEvent *evt = new MTEvent();

  TBranch *branch= tree->GetBranch("MTEvent");
  branch->SetAddress(&evt);







  for (int i=0;i<nEvent;i++){  

    if (i%100==0){cout<<i<<" \r"<<flush;}

    nevent++;
    tree.GetEntry(i);
    nchannel = evt->GetNchannel();

    nhit = 0;
    nhit_intime = 0;

    for (int j=0;j<nchannel;j++){
      
      MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);

      if (channel->GetChamberId()==5){

        nhit++;

        //if (chipid = channel->GetChipId()==tested_chip){

          xpos = channel->GetX();
          ypos = channel->GetY();

          hardid = channel->GetHardId();
      
          hxy_m2->Fill(ypos,xpos);
        
          t1 = channel->GetBcId_Abs();
          t2 = channel->GetBcId_Dif();
          t3 = channel->GetBcId_Hit();
          dt = (t2-t3);
          t = (t1 - (t2-t3));
        
          hdt->Fill(dt);
        
          if (dt>=6 && dt<=7){
          
            nhit_intime++;
            nhit_intime_tot++;
        
            hxy_m2_cut_px->Fill(ypos);
            hxy_m2_cut_py->Fill(xpos);
            hxy_m2_cut->Fill(ypos,xpos);}}}

    hnhit->Fill(nhit);
    hnhit_intime->Fill(nhit_intime);}}



  double noise =  0;
  double snr = 0;

  noise =  hdt->Integral(20,400)*2./380.;
  signal = hdt->Integral(17,18)-noise;
  snr = floor(signal/noise);

  cout<<endl;
  cout<<"Signal = "<<signal<<endl;
  cout<<"Noise = "<<noise<<endl;
  cout<<"Signal to noise ratio = "<<snr<<endl;
  cout<<"Raw eff = "<<(hdt->Integral(17,18)-noise*2)/hdt->GetEntries()*100.<<" %"<<endl;
  cout<<"SNR = "<<snr<<endl;


}

---