/Users/jacquem/NetBeansProjects/MicromegasFramework/trunk/src/analyse/Gassiplex/Shower_profile/Profile_energy.C

Go to the documentation of this file.

{gROOT->Reset();
  
  gSystem->Load("/home1/chefdevi/micromegasFrameWork/trunk/libMicro.so");
  

  bool display = 0;
  bool write_root_file = 0;


  TString fileName= "electron_2gev_4absorber.root";

  
  TString path = "/lapp_data/LC/Detecteurs/MicroMegas/data/TB2009/May/Rootfiles/dev/";
  //TString path = "/home1/chefdevi/micromegasFrameWork/trunk/src/analyse/Gassiplex/Shower_profile/";
  
  TFile *f = new TFile(path+fileName);
  cout<<"reading File : "<<path<<fileName<<endl;
  TTree *tree;
  tree  = (TTree*) f.Get("t1");cout<<"tree built\n";
  int nEvent = tree.GetEntries();
  MTEvent *evt = new MTEvent();
  TBranch *branch= tree->GetBranch("MTEvent");
  branch->SetAddress(&evt);
  cout <<"nEvent="<<nEvent<<endl;
  int NC = 672;



  /*Rootfile to save the histos*/

  if (write_root_file){
    TString fileName2= "Root_histo_files/adc_" + fileName;
    TFile * tf = new TFile(fileName2,"recreate");}



  /*Distributions of the energy (or ADC)*/

  TH1F * hadc1 = new TH1F("hadc1","",400,0,1024*5);
  TH1F * hadc2 = new TH1F("hadc2","",400,0,1024*5);
  TH1F * hadc3 = new TH1F("hadc3","",400,0,1024*5);
  TH1F * hadc4 = new TH1F("hadc4","",400,0,1024*10);

  TH1F * hadc1_clean = new TH1F("hadc1_clean","",400,0,1024*5);
  TH1F * hadc2_clean = new TH1F("hadc2_clean","",400,0,1024*5);
  TH1F * hadc3_clean = new TH1F("hadc3_clean","",400,0,1024*5);
  TH1F * hadc4_clean = new TH1F("hadc4_clean","",400,0,1024*10);

  TH1F * hadc1_no_sat = new TH1F("hadc1_no_sat","",400,0,1024*5);
  TH1F * hadc2_no_sat = new TH1F("hadc2_no_sat","",400,0,1024*5);
  TH1F * hadc3_no_sat = new TH1F("hadc3_no_sat","",400,0,1024*5);
  TH1F * hadc4_no_sat = new TH1F("hadc4_no_sat","",400,0,1024*10);

  TH1F * hpedestal = new TH1F("hpedestal","",2048,-1024,1024);
  TH1F * hmean_pedestal = new TH1F("hmean_pedestal","",2048,-1024,1024);
  
  

  /*Set Line Width*/

  hadc1->SetLineWidth(2);
  hadc2->SetLineWidth(2);
  hadc3->SetLineWidth(2);
  hadc4->SetLineWidth(2);

  hadc1_clean->SetLineWidth(2);
  hadc2_clean->SetLineWidth(2);
  hadc3_clean->SetLineWidth(2);
  hadc4_clean->SetLineWidth(2);

  hadc1_no_sat->SetLineWidth(2);
  hadc2_no_sat->SetLineWidth(2);
  hadc3_no_sat->SetLineWidth(2);
  hadc4_no_sat->SetLineWidth(2);

  /*Set Line Color*/

  hadc2->SetLineColor(2);
  hadc3->SetLineColor(4);
  hadc4->SetLineColor(8);

  hadc2_clean->SetLineColor(2);
  hadc3_clean->SetLineColor(4);
  hadc4_clean->SetLineColor(8);

  hadc2_no_sat->SetLineColor(2);
  hadc3_no_sat->SetLineColor(4);
  hadc4_no_sat->SetLineColor(8);


  /*Set Line Style*/
 
  hadc1_clean->SetLineStyle(2);
  hadc2_clean->SetLineStyle(2);
  hadc3_clean->SetLineStyle(2);
  hadc4_clean->SetLineStyle(2);

  hadc1_no_sat->SetLineStyle(2);
  hadc2_no_sat->SetLineStyle(2);
  hadc3_no_sat->SetLineStyle(2);
  hadc4_no_sat->SetLineStyle(2);


  /*Define variables*/

  int xpos1,ypos1,adc1,sum_adc1;
  int xpos2,ypos2,adc2,sum_adc2;
  int xpos3,ypos3,adc3,sum_adc3;
  int xpos4,ypos4,adc4,sum_adc4;

  int threshold = 35;
  float nhit1=0,nhit2=0,nhit3=0,nhit4=0;
  int nclean = 0;
  bool clean = 0;
  int adc_sat = 900;
  int pedestal = 0;
  int correction = 0;
  bool no_sat1 = 1,no_sat2 = 1;
  bool no_sat3 = 1,no_sat4 = 1;

  int event = 0;
  int crap_event = 0;





  /*Loop over all events*/

  for( int i=0;i<nEvent;i++){
    //for( int i=0;i<1;i++){
  //for( int i=event;i<event+1;i++){
    
    if (display){cout<<"Event # "<<i<<endl;}

    tree.GetEntry(i);
    NC = evt->GetNchannel();

    if (NC == 672){
    

    clean=0;
    no_sat1 = 1;no_sat2 = 1;no_sat3 = 1;no_sat4 = 1;
    nhit1=0,nhit2=0,nhit3=0,nhit4=0;
    sum_adc1=0;sum_adc2=0;sum_adc3=0;sum_adc4=0;
    hpedestal->Reset();


    /*Chamber beta 2.1 D*/

    for(int j=0;j<96;j++){

      MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);
      xpos1 = channel.GetX();
      ypos1 = channel.GetY();
      adc1  = channel.GetValue();
      
      if (adc1>=threshold){
        
        sum_adc1+=adc1;
        /*no saturation if all adcs in chamber are below 900*/
        if (adc1>adc_sat){no_sat1 = 0;}
        if (display){cout<<" x1="<<xpos1<<"  y1="<<ypos1<<"  adc1="<<adc1<<endl;}
        nhit1++;}}



    /*Chamber beta 2.1 E*/

    for(int j=0+96;j<96+96;j++){

      MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);
      xpos2 = channel.GetX();
      ypos2 = channel.GetY();
      adc2  = channel.GetValue();

      if (adc2>=threshold){

        sum_adc2+=adc2;
        /*no saturation if all adcs in chamber are below 900*/
        if (adc2>adc_sat){no_sat2 = 0;}
        if (display){cout<<" x2="<<xpos2<<"  y2="<<ypos2<<"  adc2="<<adc2<<endl;}
        nhit2++;}}


    /*Chamber beta 2.1 A*/

    for(int j=0+96+96;j<96+96+96;j++){

      MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);
      xpos3 = channel.GetX();
      ypos3 = channel.GetY();
      adc3  = channel.GetValue();

      if (adc3>=threshold){

        sum_adc3+=adc3;
        /*no saturation if all adcs in chamber are below 900*/
        if (adc3>adc_sat){no_sat3 = 0;}
        if (display){cout<<" x3="<<xpos3<<"  y3="<<ypos3<<"  adc3="<<adc3<<endl;}
        nhit3++;}}


    

    /*clean event if one and only one hit in each beta 2.1 chamber*/
    clean = (nhit1==1) && (nhit2==1) && (nhit3==1);
    


    /*chamber beta 2.4*/

    for(int j=287;j<672;j++){

      MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);
      xpos4 = channel.GetX();
      ypos4 = channel.GetY();
      adc4  = channel.GetValue();

      if (adc4>=threshold){
        
        sum_adc4+=adc4;

        /*no saturation if all adcs in chamber are below 900*/
        if (adc4>adc_sat){no_sat4 = 0;}
        if (display){cout<<" x4="<<xpos4<<"  y4="<<ypos4<<"  adc4="<<adc4<<endl;}
        nhit4++;}

      else{hpedestal->Fill(adc4);}}


    pedestal = hpedestal->GetMean();
    hmean_pedestal->Fill(pedestal);
    correction = nhit4 * pedestal;

    
    if (sum_adc1>0){hadc1->Fill(sum_adc1-correction);}
    if (sum_adc2>0){hadc2->Fill(sum_adc2-correction);}
    if (sum_adc3>0){hadc3->Fill(sum_adc3-correction);}
    if (sum_adc4>0){hadc4->Fill(sum_adc4-correction);}
    
    
    if (clean){
      
      if (sum_adc1>0){hadc1_clean->Fill(sum_adc1-correction);}
      if (sum_adc2>0){hadc2_clean->Fill(sum_adc2-correction);}
      if (sum_adc3>0){hadc3_clean->Fill(sum_adc3-correction);}
      if (sum_adc4>0){hadc4_clean->Fill(sum_adc4-correction);}
    
      if (no_sat1 && (sum_adc1>0)){hadc1_no_sat->Fill(sum_adc1-correction);}
      if (no_sat2 && (sum_adc2>0)){hadc2_no_sat->Fill(sum_adc2-correction);}
      if (no_sat3 && (sum_adc3>0)){hadc3_no_sat->Fill(sum_adc3-correction);}
      if (no_sat4 && (sum_adc4>0)){hadc4_no_sat->Fill(sum_adc4-correction);}

      if (display){cout<<"clean event!"<<endl;}

      nclean++;}
      

    if (display){cout<<nhit1<<" "<<nhit2<<" "<<nhit3<<" "<<nhit4<<endl;}}


    else{crap_event++;}}


      



  cout<<"ntotal = "<<nEvent<<"  nclean = "<<nclean<<"  ncrap = "<<crap_event<<endl;



  /*save histo in a rootfile*/
  if (write_root_file){
    tf->Write();
    tf->Close();}


}

---