/data3/calcul/jacquem/working_dir/Micromegas/micromegasFrameWork/src/analyse/Microroc/Calibration/Align_pedestals_asu_microroc.cpp

Go to the documentation of this file.

#include <map>
#include <fstream>
#include "TKey.h"
#include <TROOT.h>
#include <TRint.h>
#include <TKey.h>
#include <TTree.h>
#include <TFile.h>
#include <TSystem.h>
#include <iostream>
#include <TGraphErrors.h>
#include <TF1.h>
#include <TH1F.h>
#include <TCanvas.h>
#include <sstream>
#include "Log.hh"

#include "root/MTRun.hh"
#include "root/MTEvent.hh"
#include "root/MTChannel.hh"
#include "root/MTMicrorocChip.hh"

#include "MicroException.hh"


using namespace std;


int main(int argc, char**argv){

  if (argc!=6){
    
    cout<<"/*****************************************************************************************/"<<endl;
    cout<<"/* ARGUMENTS TO THE PROGRAM **************************************************************/"<<endl;
    cout<<"/* 1 - DIRECTORY OF THE 5 TEXT FILES WITH SCURVE PARAMETERS FOR THE 5 DIFFERENT OFFSETS  */"<<endl;
    cout<<"/* 2 - ASU NUMBER: 1314, 1516, 1718, 2021, 2223, 2425 etc...                             */"<<endl;
    cout<<"/* 3 - RUN NUMBER                                                                        */"<<endl;
    cout<<"/* 4 - VERBOSE (0/1)                                                                     */"<<endl;
    cout<<"/* 5 - WRITE TO FILE (0/1)                                                               */"<<endl;
    cout<<"/*****************************************************************************************/"<<endl;}

  else{
   
    
    TString dir  = argv[1];
    TString asu  = argv[2];
    TString run  = argv[3];
    int verbose = atoi(argv[4]);
    int writetofile = atoi(argv[5]);

    int nsigma = 7;

    //Scurve parameters files
    TString file0   = dir + "scurve_noise_offset0_dac0_run"  + run + "_asu" + asu + ".txt";
    TString file3   = dir + "scurve_noise_offset3_dac0_run"  + run + "_asu" + asu + ".txt";
    TString file7   = dir + "scurve_noise_offset7_dac0_run"  + run + "_asu" + asu + ".txt";
    TString file11  = dir + "scurve_noise_offset11_dac0_run" + run + "_asu" + asu + ".txt";
    TString file15  = dir + "scurve_noise_offset15_dac0_run" + run + "_asu" + asu + ".txt";

    //Noisy channel files
    TString noisy0  = dir + "noisy_noise_offset0_dac0_run"  + run + "_asu" + asu + ".txt";
    TString noisy3  = dir + "noisy_noise_offset3_dac0_run"  + run + "_asu" + asu + ".txt";
    TString noisy7  = dir + "noisy_noise_offset7_dac0_run"  + run + "_asu" + asu + ".txt";
    TString noisy11 = dir + "noisy_noise_offset11_dac0_run" + run + "_asu" + asu + ".txt";
    TString noisy15 = dir + "noisy_noise_offset15_dac0_run" + run + "_asu" + asu + ".txt";

    //Output files
    TString outfile1 = dir + "ASU" + asu + "_seuils_dac0_endpoint_align_run" + run + ".txt";
    TString outfile2 = dir + "ASU" + asu + "_gains_endpoint_align_run"  + run + ".txt";
    TString outfile3 = dir + "ASU" + asu + "_voies_bruyantes_endpoint_align_run"  + run + ".txt";

    ofstream out1(outfile1.Data());
    ofstream out2(outfile2.Data());
    ofstream out3(outfile3.Data());

    ifstream in0(file0);
    ifstream in3(file3);
    ifstream in7(file7);
    ifstream in11(file11);
    ifstream in15(file15);

    ifstream innoise0(noisy0);
    ifstream innoise3(noisy3);
    ifstream innoise7(noisy7);
    ifstream innoise11(noisy11);
    ifstream innoise15(noisy15);

 
    TGraphErrors * tg_endpoint[144][64];

    bool asu_chip[144];
    bool noisy_channel[144][64];

    float endpoint_min_15[144];
    float endpoint_max_0[144];
    float endpoint_target[144];

    TH1F * hinflex[144];
    TH1F * hwidth[144];
    TH1F * hendpoint[144];

    for (int i=0;i<144;i++)
      {
        asu_chip[i] = 0;

        endpoint_min_15[i] = 1025;
        endpoint_max_0[i] = 0;

        endpoint_target[i] = 0;

        TString name = Form("hinflex_%i",i+1);
        TString title = Form("Inflexion at offset=7 - chip %i;inflexion point (DAC)",i+1);
        hinflex[i] = new TH1F(name,title,2000,100,300);

        name = Form("hwidth_%i",i+1);
        title = Form("Width at offset=7 - chip %i;width (DAC)",i+1);
        hwidth[i] = new TH1F(name,title,2000,0,5);

        name = Form("hendpoint_%i",i+1);
        title = Form("Endpoint (%i sigma) at offset=7 - chip %i;endpoint (DAC)",nsigma,i+1);
        hendpoint[i] = new TH1F(name,title,2000,100,300);

        for (int j=0;j<64;j++)
          {
            noisy_channel[i][j] = 0;

            tg_endpoint[i][j] = new TGraphErrors();

            name = Form("tg_%i_%i",i+1,j);
            title = Form("Endpoint versus offset - chip %i - channel %i",i+1,j);            

            tg_endpoint[i][j]->SetName(name);
            tg_endpoint[i][j]->SetTitle(title);

            tg_endpoint[i][j]->SetMarkerStyle(20);
            tg_endpoint[i][j]->SetMarkerColor(4);
          }
      }





    int N = 0;
    int offset = 0;
    int chip = 0;
    int channel = 0;
    float inflex = 0;
    float width = 0;
    float endpoint = 0;
    float inflex_err = 0;
    float width_err = 0;
    float endpoint_err = 0;


    

    /*Read noisy channels*/
    cout<<endl;
    cout<<"READ NOISY CHANNELS FILES"<<endl;
    cout<<"Read file "<<noisy0<<endl; 
    while(innoise0>>chip>>channel>>inflex>>inflex_err>>width>>width_err) {noisy_channel[chip-1][channel] = 1;}
    cout<<"Read file "<<noisy3<<endl; 
    while(innoise3>>chip>>channel>>inflex>>inflex_err>>width>>width_err) {noisy_channel[chip-1][channel] = 1;}
    cout<<"Read file "<<noisy7<<endl; 
    while(innoise7>>chip>>channel>>inflex>>inflex_err>>width>>width_err) {noisy_channel[chip-1][channel] = 1;}
    cout<<"Read file "<<noisy11<<endl; 
    while(innoise11>>chip>>channel>>inflex>>inflex_err>>width>>width_err){noisy_channel[chip-1][channel] = 1;}
    cout<<"Read file "<<noisy15<<endl; 
    while(innoise15>>chip>>channel>>inflex>>inflex_err>>width>>width_err){noisy_channel[chip-1][channel] = 1;}


    //Read Scurve files - Offset = 0
    cout<<endl;
    cout<<"READ SCURVE PARAMETERS FILES"<<endl;
    cout<<"Read file "<<file0<<endl;
    offset = 0;
    while(in0>>chip>>channel>>inflex>>inflex_err>>width>>width_err)
      {
        if (!noisy_channel[chip-1][channel])
          {
            asu_chip[chip-1] = 1;

            endpoint = inflex + nsigma*width;
            endpoint_err = sqrt(inflex_err + nsigma*nsigma*width_err);

            N = tg_endpoint[chip-1][channel]->GetN();

            tg_endpoint[chip-1][channel]->SetPoint(N,offset,endpoint);
            tg_endpoint[chip-1][channel]->SetPointError(N,0,endpoint_err);

            if (endpoint>endpoint_max_0[chip-1]){endpoint_max_0[chip-1] = endpoint;}
          }
      }
 

    //Read Scurve files - Offset = 3
    cout<<"Read file "<<file3<<endl;
    offset = 3;
    while(in3>>chip>>channel>>inflex>>inflex_err>>width>>width_err)
      {
        if (!noisy_channel[chip-1][channel])
          {
            endpoint = inflex + nsigma*width;
            endpoint_err = sqrt(inflex_err + nsigma*nsigma*width_err);

            N = tg_endpoint[chip-1][channel]->GetN();

            tg_endpoint[chip-1][channel]->SetPoint(N,offset,endpoint);
            tg_endpoint[chip-1][channel]->SetPointError(N,0,endpoint_err);
          }
      }


    //Read Scurve files - Offset = 7
    cout<<"Read file "<<file7<<endl;
    offset = 7;
    while(in7>>chip>>channel>>inflex>>inflex_err>>width>>width_err)
      {
        if (!noisy_channel[chip-1][channel])
          {
            endpoint = inflex + nsigma*width;
            endpoint_err = sqrt(inflex_err + nsigma*nsigma*width_err);

            hinflex[chip-1]->Fill(inflex);
            hwidth[chip-1]->Fill(width);
            hendpoint[chip-1]->Fill(inflex+nsigma*width);

            N = tg_endpoint[chip-1][channel]->GetN();

            tg_endpoint[chip-1][channel]->SetPoint(N,offset,endpoint);
            tg_endpoint[chip-1][channel]->SetPointError(N,0,endpoint_err);
          }
      }


    //Read Scurve files - Offset = 11
    cout<<"Read file "<<file11<<endl;
    offset = 11;
    while(in11>>chip>>channel>>inflex>>inflex_err>>width>>width_err)
      {
        if (!noisy_channel[chip-1][channel])
          {
            endpoint = inflex + nsigma*width;
            endpoint_err = sqrt(inflex_err + nsigma*nsigma*width_err);

            N = tg_endpoint[chip-1][channel]->GetN();

            tg_endpoint[chip-1][channel]->SetPoint(N,offset,endpoint);
            tg_endpoint[chip-1][channel]->SetPointError(N,0,endpoint_err);
          }
      }


    //Read Scurve files - Offset = 15
    cout<<"Read file "<<file15<<endl;
    offset = 15;
    while(in15>>chip>>channel>>inflex>>inflex_err>>width>>width_err)
      {
        if (!noisy_channel[chip-1][channel])
          {
            endpoint = inflex + nsigma*width;
            endpoint_err = sqrt(inflex_err + nsigma*nsigma*width_err);

            N = tg_endpoint[chip-1][channel]->GetN();

            tg_endpoint[chip-1][channel]->SetPoint(N,offset,endpoint);
            tg_endpoint[chip-1][channel]->SetPointError(N,0,endpoint_err);

            if (endpoint<endpoint_min_15[chip-1]){endpoint_min_15[chip-1] = endpoint;}
          }
      }




    //Determine endpoint onto which to align Scurves
    for (int i=0;i<144;i++)
      {
        if (asu_chip[i])
          {
            endpoint_target[i] = 0.5 * (endpoint_max_0[i] + endpoint_min_15[i]);

            bool alignment_possible = 0;
            if (endpoint_max_0[i]<=endpoint_min_15[i])
              {
                alignment_possible = 1;
              }

            if(verbose == 1)
              {
                cout<<endl;
                cout<<"Chip "<<i+1<<endl;
                cout<<"  Endpoint max at offset 0 = "<<endpoint_max_0[i]<<endl;
                cout<<"  Endpoint min at offset 15 = "<<endpoint_min_15[i]<<endl;

                if (alignment_possible)
                  {
                    cout<<"   -> Endpoint alignment possible ******************"<<endl;
                    cout<<"      ->  Will be aligned on "<<floor(endpoint_target[i]+1)<<endl;
                  }
                else
                  {
                    cout<<"   -> Endpoint alignment NOT possible"<<endl;
                  }
              }
          }
      }



    //Calculate offset map and thresholds

    float a = 0;
    float b = 0;
    float offset_ped = 0;
    float offset_end = 0;

    TF1 * f = new TF1("f","pol1",0,14);

    for (int i=0;i<144;i++)
      {
        if (hinflex[i]->GetEntries()!=0)
          {
            //Write thresholds
            out1<<i+1<<" "<<floor(endpoint_target[i] + 1)<<endl;

            //Write offsets and noisy channels
            for (int j=0;j<64;j++)
              {
                if (noisy_channel[i][j])
                  {
                    out2<<i+1<<" "<<j<<" "<<0<<endl;
                    out3<<i+1<<" "<<j<<endl;
                  }
                else
                  {
                    //Endpoint versus offset
                    tg_endpoint[i][j]->Fit(f,"q");

                    a = f->GetParameter(1);
                    b = f->GetParameter(0);

                    if (a!=0)
                      {
                        offset_end = floor((endpoint_target[i] - b)/a + 0.5);
                      }
                    else
                      {
                        offset_end =  7;
                      }

                    if (offset_end<0)
                      {
                        offset_end = 0;
                      }
                    if (offset_end>15)
                      {
                        offset_end = 15;
                      }

                    out2<<i+1<<" "<<j<<" "<<offset_end<<endl;
                  }
              }
          }
      }
 

    //Write channel TGraph and chip histos to ROOT file
    if (writetofile == 1)
      {
        TString rootfilename = dir + "pedestal_alignment.root";

        cout<<endl;
        cout<<"WRITE CURVES AND HISTOS TO ROOT FILE : "<<rootfilename<<endl;

        TFile * tf = new TFile(rootfilename,"recreate");

        for (int i=0;i<144;i++)
          {
            if (hinflex[i]->GetEntries() != 0)
              {
                hinflex[i]->Write();
                hwidth[i]->Write();
                hendpoint[i]->Write();

                for (int j=0;j<64;j++)
                  {
                    tg_endpoint[i][j]->Write();
                  }
              }
          }
        tf->Close();
      }



  }
}

---