Changeset 208 in JEM-EUSO


Ignore:
Timestamp:
Sep 30, 2013, 10:04:31 AM (11 years ago)
Author:
moretto
Message:

addition of default gain tables for ASIC C & D needed to create the default gain-org.txt gain table

Location:
equalization_gain/trunk
Files:
2 added
5 edited

Legend:

Unmodified
Added
Removed
  • equalization_gain/trunk/README_camille

    r191 r208  
    108108% .q
    109109
     1108bis) You can also run the with root the App PlotScurvesAllMain as follow:
     111
     112% ../../src/PlotScurvesAllMain
     113
     114It will also produce the .data file described at step 9). You have to go directly to step 10) then.
     115
     116However you will have to "make clean" & "make" after changing fitting parameters etc.
     117
    1101189) Now we chose the header parameters we can make the new gain table. Run:
    111119
     
    117125Do the same instructions for an other ASIC (D one for example).
    118126
    119 Go back in the previous folder to run the makeGainTable.sh bash script. But first, modify parameters of it as name of file and created .txt. You also have to make modification depending on whether you use type C or D ASICs :
     12710) Go back in the previous folder to run the makeGainTable.sh bash script. But first, modify parameters of it as name of file and created .txt. You also have to make modification depending on whether you use type C or D ASICs :
    120128
    121129% cd ../
  • equalization_gain/trunk/src/PlotScurvesAll-diffFit.cc

    r203 r208  
    1818#include "PlotScurvesAll-diffFit.h"
    1919
    20 //#if defined(__CINT__) && !defined(__MAKECINT__) 
     20//#if defined(__CINT__) && !defined(__MAKECINT__)
    2121//#include "PlotScurveAll-diffFit.h"
    22 //#else 
     22//#else
    2323//#include "PlotScurvesAll-diffFit.h"
    2424//#endif
     
    4040void Help()
    4141{
    42 
    43   cout << "This is the Help of the PlotScurveAll-diffFit" <<endl;
    44   cout << " Version 1.0 frol September 24 2013 " <<endl;
    45   cout << " Main author Hiroko Miyamoto " << endl;
    46   cout << " ==> Please in ROOT interpreter, issue the command .L PlotScurvesAll-diffFit.cc " <<endl;
    47   cout << " ==> Then play with PlotScurve() function " << endl;
    48   cout << " ==> If you need more information, please type Menu()" << endl;
    49    
     42   
     43    cout << "This is the Help of the PlotScurveAll-diffFit" <<endl;
     44    cout << " Version 1.0 frol September 24 2013 " <<endl;
     45    cout << " Main author Hiroko Miyamoto " << endl;
     46    cout << " ==> Please in ROOT interpreter, issue the command .L PlotScurvesAll-diffFit.cc " <<endl;
     47    cout << " ==> Then play with PlotScurve() function " << endl;
     48    cout << " ==> If you need more information, please type Menu()" << endl;
     49    
    5050}
    5151//--------------------------------------------------------------------------------------------------------
     
    5454void Menu()
    5555{
    56  cout << "This is the Detailed menu information for running the PlotScurveAll-diffFit" <<endl;
    57 
    58  cout << " void PlotScurve( " << endl;
    59  cout << " \t " <<  " Int_t    kDrawCh=31, " << endl;
    60  cout << " \t " <<  " TString  fname=kFname, " << endl;
    61  cout << " \t " <<  " TString  fname_ped=kFnamePed, //ASIC_D ped_file " << endl;
    62  cout << " \t " <<  " Double_t dac_half_spe_fix = -1," << endl;
    63  cout << " \t " <<  " TString  fname_gain=\"gain-org.txt\"," << endl;
    64  cout << " \t " <<  " Int_t    fit_min=kFitMin," << endl;
    65  cout << " \t " <<  " Int_t    fit_max=kFitMax," << endl;
    66  cout << " \t " <<  " Int_t    fit_min_scurve=kFitMinS," << endl;
    67  cout << " \t " <<  " Int_t    fit_max_scurve=kFitMaxS," << endl;
    68  cout << " \t " <<  " Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D " << endl;
    69  cout << " \t " <<  " Int_t    kChRef=kRefCh,//ASIC_C " << endl;
    70  cout << " \t " <<  " Bool_t   checkAve=0 //0:use ref pix, 1:check average " << endl;
    71  cout << " ) " << endl;         
    72  cout << "Set the parameters on the top of the program depending on your analysis:  " << endl;
    73  cout << "General policy: parameters starting with a \"k\" at the beginning is fixed gloval parameter and usually set in the header file, never change through the program. " << endl;
    74  cout << " : " << endl;
    75  cout << " - kDrawCh: channel just you want to check on display. doens't effect the analysis. " << endl;
    76  
    77  cout << " - fname, fname_ped (kFname, kFnamePed): \"YOUR DATA FILE\" without .data"  << endl;
    78  cout << " See examples in the header file. " << endl;
    79  cout << " - dac_half_spe_fix: set DAC value only when you want to adjust the gain to particular pulse height (DAC value) instead of either average of the pulse height of reference pixel. Usually just leave it -1. " << endl;
    80  cout << " - fit_min, fit_max (kFitMin, kFitMax): leave them as default, it is not currently used for the analysis."  << endl;
    81  cout << " - fit_min_scurve, fit_max_scurve (kFitMinS, kFitMaxS): define the range of s-curve fitting. search proper values by yourself with running root program a few times."  << endl;
    82  cout << " - kBoard: don't forget to change 0(C) or 1(D) depending on which ASIC you're analysing "  << endl;
    83  cout << " - kChRef (kRefCh): a pixel whose gain remains 64 thorough the gain adjustment. " << endl;
    84  cout << " - checkAve: for the first run, set this parameter 1, and find a pixel which remains 64. In the case of parameter 1, all the gains are adjusted to the average pulse height. " << endl;
    85    cout <<  " Once you choose the reference pixel, set the parameter 0, as well as set kRefCh the pixel number. (better to keep in the header file for the record of each PMTs.) " << endl;
    86 
    87    cout << " Actual configuration " << endl;
    88    cout << " \t " <<  "TString kFname = " << kFname << endl;
    89    cout << " \t " <<  "TString  kFnamePed = " << kFnamePed  << endl;
    90    cout << " \t " <<  "Int_t kFitMin = " << kFitMin << endl;
    91    cout << " \t " <<  "Int_t kFitMax = " << kFitMax << endl;
    92 
     56    cout << "This is the Detailed menu information for running the PlotScurveAll-diffFit" <<endl;
     57   
     58    cout << " void PlotScurve( " << endl;
     59    cout << " \t " <<  " Int_t    kDrawCh=31, " << endl;
     60    cout << " \t " <<  " TString  fname=kFname, " << endl;
     61    cout << " \t " <<  " TString  fname_ped=kFnamePed, //ASIC_D ped_file " << endl;
     62    cout << " \t " <<  " Double_t dac_half_spe_fix = -1," << endl;
     63    cout << " \t " <<  " TString  fname_gain=\"gain-org.txt\"," << endl;
     64    cout << " \t " <<  " Int_t    fit_min=kFitMin," << endl;
     65    cout << " \t " <<  " Int_t    fit_max=kFitMax," << endl;
     66    cout << " \t " <<  " Int_t    fit_min_scurve=kFitMinS," << endl;
     67    cout << " \t " <<  " Int_t    fit_max_scurve=kFitMaxS," << endl;
     68    cout << " \t " <<  " Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D " << endl;
     69    cout << " \t " <<  " Int_t    kChRef=kRefCh,//ASIC_C " << endl;
     70    cout << " \t " <<  " Bool_t   checkAve=0 //0:use ref pix, 1:check average " << endl;
     71    cout << " ) " << endl;
     72    cout << "Set the parameters on the top of the program depending on your analysis:  " << endl;
     73    cout << "General policy: parameters starting with a \"k\" at the beginning is fixed gloval parameter and usually set in the header file, never change through the program. " << endl;
     74    cout << " : " << endl;
     75    cout << " - kDrawCh: channel just you want to check on display. doens't effect the analysis. " << endl;
     76    
     77    cout << " - fname, fname_ped (kFname, kFnamePed): \"YOUR DATA FILE\" without .data"  << endl;
     78    cout << " See examples in the header file. " << endl;
     79    cout << " - dac_half_spe_fix: set DAC value only when you want to adjust the gain to particular pulse height (DAC value) instead of either average of the pulse height of reference pixel. Usually just leave it -1. " << endl;
     80    cout << " - fit_min, fit_max (kFitMin, kFitMax): leave them as default, it is not currently used for the analysis."  << endl;
     81    cout << " - fit_min_scurve, fit_max_scurve (kFitMinS, kFitMaxS): define the range of s-curve fitting. search proper values by yourself with running root program a few times."  << endl;
     82    cout << " - kBoard: don't forget to change 0(C) or 1(D) depending on which ASIC you're analysing "  << endl;
     83    cout << " - kChRef (kRefCh): a pixel whose gain remains 64 thorough the gain adjustment. " << endl;
     84    cout << " - checkAve: for the first run, set this parameter 1, and find a pixel which remains 64. In the case of parameter 1, all the gains are adjusted to the average pulse height. " << endl;
     85    cout <<  " Once you choose the reference pixel, set the parameter 0, as well as set kRefCh the pixel number. (better to keep in the header file for the record of each PMTs.) " << endl;
     86   
     87    cout << " Actual configuration " << endl;
     88    cout << " \t " <<  "TString kFname = " << kFname << endl;
     89    cout << " \t " <<  "TString  kFnamePed = " << kFnamePed  << endl;
     90    cout << " \t " <<  "Int_t kFitMin = " << kFitMin << endl;
     91    cout << " \t " <<  "Int_t kFitMax = " << kFitMax << endl;
     92   
    9393}
    9494//--------------------------------------------------------------------------------------------------------
     
    9797void InitRoot()
    9898{
    99 
    100   gStyle->SetTitleFillColor(kWhite);
    101   gStyle->SetTitleColor(kBlue,"xyz");
    102   gStyle->SetTitleSize(0.03,"xyz");
    103   gStyle->SetLabelSize(0.03,"xyz");
    104   gStyle->SetLabelColor(kBlue,"xyz");
    105   gStyle->SetAxisColor(kMagenta);
    106   gStyle->SetPalette(1,0);
    107   gStyle->SetOptFit(101);
    108   gStyle->SetOptStat(0);
    109  
    110   //  gStyle->SetOptStat("emr");
    111   gStyle->SetOptStat("");
     99   
     100    gStyle->SetTitleFillColor(kWhite);
     101    gStyle->SetTitleColor(kBlue,"xyz");
     102    gStyle->SetTitleSize(0.03,"xyz");
     103    gStyle->SetLabelSize(0.03,"xyz");
     104    gStyle->SetLabelColor(kBlue,"xyz");
     105    gStyle->SetAxisColor(kMagenta);
     106    gStyle->SetPalette(1,0);
     107    gStyle->SetOptFit(101);
     108    gStyle->SetOptStat(0);
     109   
     110    //  gStyle->SetOptStat("emr");
     111    gStyle->SetOptStat("");
    112112}
    113113//-------------------------------------------------------------------
     
    116116//--------------------------------------------------------------------------------------------------------
    117117/*
    118 
    119 The declaration of this function is now put inside the h file with the default values to be used
    120 
    121 void PlotScurve(
    122                 Int_t    kDrawCh=31,
    123                 TString  fname=kFname,
    124                 TString  fname_ped=kFnamePed, //ASIC_D ped_file
    125                 Double_t dac_half_spe_fix = -1,
    126                 TString  fname_gain="gain-org.txt",
    127                 Int_t    fit_min=kFitMin,
    128                 Int_t    fit_max=kFitMax,
    129                 Int_t    fit_min_scurve=kFitMinS,
    130                 Int_t    fit_max_scurve=kFitMaxS,
    131                 Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D
    132                 Int_t    kChRef=kRefCh,//ASIC_C
    133                 Bool_t   checkAve=0 //0:use ref pix, 1:check average
    134                 )
    135 */
     118 
     119 The declaration of this function is now put inside the h file with the default values to be used
     120 
     121 void PlotScurve(
     122 Int_t    kDrawCh=31,
     123 TString  fname=kFname,
     124 TString  fname_ped=kFnamePed, //ASIC_D ped_file
     125 Double_t dac_half_spe_fix = -1,
     126 TString  fname_gain="gain-org.txt",
     127 Int_t    fit_min=kFitMin,
     128 Int_t    fit_max=kFitMax,
     129 Int_t    fit_min_scurve=kFitMinS,
     130 Int_t    fit_max_scurve=kFitMaxS,
     131 Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D
     132 Int_t    kChRef=kRefCh,//ASIC_C
     133 Bool_t   checkAve=0 //0:use ref pix, 1:check average
     134 )
     135 */
    136136
    137137// Definition of the function  PlotScurve
    138138void PlotScurve(
    139                 Int_t    kDrawCh,
    140                 TString  fname,
    141                 TString  fname_ped, //ASIC_D ped_file
    142                 Double_t dac_half_spe_fix,
    143                 TString  fname_gain,
    144                 Int_t    fit_min,
    145                 Int_t    fit_max,
    146                 Int_t    fit_min_scurve,
    147                 Int_t    fit_max_scurve,
    148                 Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
    149                 Int_t    kChRef,//ASIC_C
    150                 Bool_t   checkAve //0:use ref pix, 1:check average
    151                 )
     139                Int_t    kDrawCh,
     140                TString  fname,
     141                TString  fname_ped, //ASIC_D ped_file
     142                Double_t dac_half_spe_fix,
     143                TString  fname_gain,
     144                Int_t    fit_min,
     145                Int_t    fit_max,
     146                Int_t    fit_min_scurve,
     147                Int_t    fit_max_scurve,
     148                Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
     149                Int_t    kChRef,//ASIC_C
     150                Bool_t   checkAve //0:use ref pix, 1:check average
     151                )
    152152
    153153
    154154{
    155                 //NOTE   
    156                 //Double_t dac_half_spe_fix = 177; //for EC_ASICb N101(asic N158).
    157                 //check DAC value corresponds to voltage to 160fC input
    158                 //(depends on DAC linearity)
    159                 //Double_t dac_half_spe_fix = -1 //set gain to DAC average (default)
    160 
    161   /******************* Readout data file *********************/
    162  
    163   /*********************** SIGNAL RUN ********************/
    164   TString *input             = new TString[kNch];
    165   TString *input_ped         = new TString[kNch];
    166   Int_t    ch;
    167   Float_t  dac_half_sum      = 0;
    168   Float_t  dac_ped_sum       = 0;
    169   Float_t  dac_spe_sum       = 0;
    170   Float_t  dac_sfit_sum      = 0;
    171   Float_t  diff_sum          = 0;
    172   Float_t  mean_sum          = 0;
    173   Float_t  maxBinSum         = 0;
    174   Double_t dac_half[kNch];
    175   Double_t dac_ped[kNch];
    176   Double_t dac_spe[kNch];
    177   Double_t dac_sfit[kNch];
    178   Double_t mean_diff[kNch];
    179   Double_t gain_org[kNch];
    180   Double_t sfit_chiS[kNch];
    181   Double_t flag_fit[kNch];
    182 
    183 
    184   InitRoot();
    185 
    186   TCanvas *c1                = new TCanvas(kCanvasName1,kCanvasTitle1,200,20,1100,600); //smoothed S-curves
    187   c1->Divide(11,6);
    188   TCanvas *c2                = new TCanvas(kCanvasName2,kCanvasTitle2,300,50,1100,600); //differentiated histograms of smoothed s-curves
    189   c2->Divide(11,6);
    190   TCanvas *c5                = new TCanvas(kCanvasName5,kCanvasTitle5,400,70,1100,600);
    191   c5->Divide(11,6);
    192   TCanvas *c3                = new TCanvas(kCanvasName3,kCanvasTitle3,500,100,1200,300);
    193   c3->Divide(4,1);
    194   /*TCanvas *c6                = new TCanvas("c6","c6",450,150,1000,500);
    195   c6->Divide(4,2);
    196   TCanvas *c7                = new TCanvas("c7","c7",470,200,1000,500);
    197   c7->Divide(4,2);*/
    198   TCanvas *c4                = new TCanvas(kCanvasName4,kCanvasTitle4,900,400,300,300);
    199 
    200   TString  name_hist[kNch];
    201   TString  name_hist_ped[kNch];
    202   TString  name_hist_diff[kNch];
    203   TString  name_hist_diff_fit[kNch];
    204   TString  name_fit_diff[kNch];
    205   TH1D    *hist_scurve       = new TH1D[kNch]; // smoothed S-curves
    206   TH1D    *hist_scurve_ped   = new TH1D[kNch]; // S-curves for pedestal runs
    207   TH1D    *hist_diff         = new TH1D[kNch]; // differentiated histograms of smoothed s-curves
    208   TH1D    *hist_diff_fit     = new TH1D[kNch]; // differentiated histograms of fitted s-curves
    209   TF1     *fit_sdiff         = new TF1[kNch];  // function fitted on differentiated histograms of fitted s-curves
    210 
    211 
    212   // Open the ASIC Gain file
    213   //-------------------------
    214 
    215   TString input_gain(fname_gain);
    216   ifstream in_gain(input_gain.Data());
    217   if(!in_gain.good()){
    218     cerr << "File " << input_gain << " open failed!" << endl;
    219     return;
    220   }
    221  
    222  
    223 
    224   // first loop on channels to fill histonames
    225   //----------------------
    226   for(ch=0;ch<kNch;ch++){
    227     dac_half[ch]=0;
    228     dac_ped[ch]=0;
    229     dac_spe[ch]=0;
    230     dac_sfit[ch]=0;
    231     gain_org[ch]=0;
    232     mean_diff[ch]=0;
    233     sfit_chiS[ch]=0;
    234     flag_fit[ch]=0;
    235     name_hist[ch]           ="hist";
    236     name_hist[ch]          +=ch;
    237     name_hist_ped[ch]       ="hist_ped";
    238     name_hist_ped[ch]      +=ch;
    239     name_hist_diff[ch]      ="hist_diff";
    240     name_hist_diff[ch]     +=ch;
    241     name_hist_diff_fit[ch]  ="hist_diff_fit";
    242     name_hist_diff_fit[ch] +=ch;
    243     name_fit_diff[ch]       ="fit_diff";
    244     name_fit_diff[ch]      +=ch;
    245   }
    246  
    247 
    248   // read gain
    249   Double_t da;
    250   Int_t count=0;
    251   ch=0;
    252   while(in_gain >> da){
    253     if(!kBoard && (count>=128 && count<192)){
    254       gain_org[ch]=da;
    255       cerr << "!kBoard, count, gain_org[" << ch << "] "
    256            << !kBoard << " " << count << " " << gain_org[ch] << endl;
    257       ch++;
    258     }
    259     if(kBoard && (count>=192 && count<256)){
    260       gain_org[ch]=da;
    261       cerr << "kBoard, count, gain_org[" << ch << "] "
    262            << kBoard << " " << count << " " << gain_org[ch] << endl;
    263       ch++;
    264     }
    265     count++;
    266   }// end if first loop on all channels     
    267   in_gain.close();
    268   in_gain.clear();
    269  
    270   //----
    271  
    272   // Open the pedestal file and the S-Curve files for each channel and test if the file can be open
    273 
    274   Int_t trigger=0, trigger_fit=0, trigger_mean=0, trigger_sfit=0;
    275   // Big loop on all the channels : second loop
    276   //--------------------------------------------
    277   // build the input data file name for each channel fname-chXX.dat
    278   for(ch=0;ch<kNch;ch++){
     155    //NOTE
     156    //Double_t dac_half_spe_fix = 177; //for EC_ASICb N101(asic N158).
     157    //check DAC value corresponds to voltage to 160fC input
     158    //(depends on DAC linearity)
     159    //Double_t dac_half_spe_fix = -1 //set gain to DAC average (default)
     160   
     161    /******************* Readout data file *********************/
     162   
     163    /*********************** SIGNAL RUN ********************/
     164    TString *input             = new TString[kNch];
     165    TString *input_ped         = new TString[kNch];
     166    Int_t    ch;
     167    Float_t  dac_half_sum      = 0;
     168    Float_t  dac_ped_sum       = 0;
     169    Float_t  dac_spe_sum       = 0;
     170    Float_t  dac_sfit_sum      = 0;
     171    Float_t  diff_sum          = 0;
     172    Float_t  mean_sum          = 0;
     173    Float_t  maxBinSum         = 0;
     174    Double_t dac_half[kNch];
     175    Double_t dac_ped[kNch];
     176    Double_t dac_spe[kNch];
     177    Double_t dac_sfit[kNch];
     178    Double_t mean_diff[kNch];
     179    Double_t gain_org[kNch];
     180    Double_t sfit_chiS[kNch];
     181    Double_t flag_fit[kNch];
     182   
     183   
     184    InitRoot();
     185   
     186    TCanvas *c1                = new TCanvas(kCanvasName1,kCanvasTitle1,200,20,1100,600); //smoothed S-curves
     187    c1->Divide(11,6);
     188    TCanvas *c2                = new TCanvas(kCanvasName2,kCanvasTitle2,300,50,1100,600); //differentiated histograms of smoothed s-curves
     189    c2->Divide(11,6);
     190    TCanvas *c5                = new TCanvas(kCanvasName5,kCanvasTitle5,400,70,1100,600);
     191    c5->Divide(11,6);
     192    TCanvas *c3                = new TCanvas(kCanvasName3,kCanvasTitle3,500,100,1200,300);
     193    c3->Divide(4,1);
     194    /*TCanvas *c6                = new TCanvas("c6","c6",450,150,1000,500);
     195     c6->Divide(4,2);
     196     TCanvas *c7                = new TCanvas("c7","c7",470,200,1000,500);
     197     c7->Divide(4,2);*/
     198    TCanvas *c4                = new TCanvas(kCanvasName4,kCanvasTitle4,900,400,300,300);
     199   
     200    TString  name_hist[kNch];
     201    TString  name_hist_ped[kNch];
     202    TString  name_hist_diff[kNch];
     203    TString  name_hist_diff_fit[kNch];
     204    TString  name_fit_diff[kNch];
     205    TH1D    *hist_scurve       = new TH1D[kNch]; // smoothed S-curves
     206    TH1D    *hist_scurve_ped   = new TH1D[kNch]; // S-curves for pedestal runs
     207    TH1D    *hist_diff         = new TH1D[kNch]; // differentiated histograms of smoothed s-curves
     208    TH1D    *hist_diff_fit     = new TH1D[kNch]; // differentiated histograms of fitted s-curves
     209    TF1     *fit_sdiff         = new TF1[kNch];  // function fitted on differentiated histograms of fitted s-curves
     210   
     211   
     212    // Open the ASIC Gain file
     213    //-------------------------
     214   
     215    TString input_gain(fname_gain);
     216    ifstream in_gain(input_gain.Data());
     217    if(!in_gain.good()){
     218        cerr << "File " << input_gain << " open failed!" << endl;
     219        return;
     220    }
     221   
     222   
     223   
     224    // first loop on channels to fill histonames
     225    //----------------------
     226    for(ch=0;ch<kNch;ch++){
     227        dac_half[ch]=0;
     228        dac_ped[ch]=0;
     229        dac_spe[ch]=0;
     230        dac_sfit[ch]=0;
     231        gain_org[ch]=0;
     232        mean_diff[ch]=0;
     233        sfit_chiS[ch]=0;
     234        flag_fit[ch]=0;
     235        name_hist[ch]           ="hist";
     236        name_hist[ch]          +=ch;
     237        name_hist_ped[ch]       ="hist_ped";
     238        name_hist_ped[ch]      +=ch;
     239        name_hist_diff[ch]      ="hist_diff";
     240        name_hist_diff[ch]     +=ch;
     241        name_hist_diff_fit[ch]  ="hist_diff_fit";
     242        name_hist_diff_fit[ch] +=ch;
     243        name_fit_diff[ch]       ="fit_diff";
     244        name_fit_diff[ch]      +=ch;
     245    }
     246   
     247   
     248    // read gain
     249    Double_t da;
     250    Int_t count=0;
     251    ch=0;
     252    while(in_gain >> da){
     253        if(!kBoard && (count>=128 && count<192)){
     254            gain_org[ch]=da;
     255            cerr << "!kBoard, count, gain_org[" << ch << "] "
     256            << !kBoard << " " << count << " " << gain_org[ch] << endl;
     257            ch++;
     258        }
     259        if(kBoard && (count>=192 && count<256)){
     260            gain_org[ch]=da;
     261            cerr << "kBoard, count, gain_org[" << ch << "] "
     262            << kBoard << " " << count << " " << gain_org[ch] << endl;
     263            ch++;
     264        }
     265        count++;
     266    }// end if first loop on all channels
     267    in_gain.close();
     268    in_gain.clear();
     269   
     270    //----
     271   
     272    // Open the pedestal file and the S-Curve files for each channel and test if the file can be open
     273   
     274    Int_t trigger=0, trigger_fit=0, trigger_mean=0, trigger_sfit=0;
     275    // Big loop on all the channels : second loop
     276    //--------------------------------------------
     277    // build the input data file name for each channel fname-chXX.dat
     278    for(ch=0;ch<kNch;ch++){
     279        input[ch] += fname;
     280        input[ch] += "-ch";
     281        input[ch] += ch;
     282        input[ch] += ".dat";
     283        ifstream in(input[ch].Data());
     284        if(!in.good()){
     285            cerr << "File " << input[ch] << " open failed!" << endl;
     286            return;
     287        }
     288       
     289        // build the input data file name for each channel fname_ped-chXX.dat
     290        input_ped[ch] += fname_ped;
     291        input_ped[ch] += "-ch";
     292        input_ped[ch] += ch;
     293        input_ped[ch] += ".dat";
     294        ifstream in_ped(input_ped[ch].Data());
     295        if(!in_ped.good()){
     296            cerr << "File " << input_ped[ch] << " open failed!" << endl;
     297            return;
     298        }
     299       
     300       
     301       
     302        //10DAC~10mV
     303        // determine the DAC range in dac_min, dac_max, step_dac from file input_ped
     304        // and count the number of events
     305        Double_t da, dat, dac_min, dac_max, step_dac;
     306        count = 0;
     307        while (in >> da >> dat ){
     308            if(count==0) dac_min  = da;
     309            if(count==1) step_dac = da - dac_min;
     310            dac_max=da;
     311            count++;
     312        }
     313        in.close();
     314        in.clear();
     315       
     316        /*
     317         if(ch==43 || ch==47){//for T1 ASIC_C
     318         cerr << "CHECK " << ch << " " << fit_min_scurve << endl;
     319         Int_t fit_min_scurve=140;
     320         } else {
     321         Int_t fit_min_scurve=kFitMinS;
     322         }
     323         */
     324       
     325        // correct the scale for DAC range
     326        if(kXmin<dac_min)kXmin=dac_min;
     327        const Int_t kNdata    = count;
     328        count = (fit_max_scurve-fit_min_scurve)/step_dac;
     329        const Int_t kNdataFit = count;
     330        if(ch==0){
     331            cerr << "kNdata, dac_min, dac_max, step_dac: "
     332            << kNdata << " " << dac_min << " " << dac_max << " "
     333            << step_dac << endl;
     334        }
     335       
     336        Double_t dac_min_ped, dac_max_ped, step_dac_ped;
     337        count = 0;
     338        while (in_ped >> da >> dat ){
     339            if(count==0) dac_min_ped  = da;
     340            if(count==1) step_dac_ped = da - dac_min_ped;
     341            dac_max_ped=da;
     342            count++;
     343        }
     344       
     345        in_ped.close();
     346        in_ped.clear();
     347       
     348       
     349        const Int_t kNdataPed = count;
     350        if(ch==0){
     351            cerr << "kNdataPed, dac_min_ped, dac_max_ped, step_dac_ped: "
     352            << kNdataPed << " " << dac_min_ped << " " << dac_max_ped
     353            << " " << step_dac_ped << endl;
     354        }
     355       
     356        // book histograms
     357        // SDC: 27/09/2013 hist_scurve is a pointer on an array of THD, thus the syntax is correct (no new operator)
     358       
     359       
     360        hist_scurve[ch]      = TH1D(name_hist[ch],name_hist[ch],kNdata-1,
     361                                    dac_min,dac_max);
     362        hist_scurve_ped[ch]  = TH1D(name_hist_ped[ch],name_hist_ped[ch],
     363                                    kNdataPed-1,dac_min_ped,dac_max_ped);
     364        hist_diff[ch]        = TH1D(name_hist_diff[ch],name_hist_diff[ch],
     365                                    kNdata-2,dac_min,dac_max);
     366        hist_diff_fit[ch]    = TH1D(name_hist_diff_fit[ch],name_hist_diff_fit[ch],
     367                                    kNdataFit,fit_min_scurve,fit_max_scurve);
     368       
     369       
     370       
     371       
     372       
     373        // read the data for that channel for normal data
     374        in.open(input[ch].Data());
     375        Double_t       data[kNdata], dac[kNdata], dac_fit[kNdataFit], data_fit[kNdataFit];
     376        while (in >> da >> dat ){
     377            hist_scurve[ch].Fill(da,dat); // fill the histograms of the s-curve
     378        }
     379        in.close();
     380        in.clear();
     381       
     382        // read the data for pedestal data
     383        in_ped.open(input_ped[ch].Data());
     384        while (in_ped >> da >> dat ){
     385            hist_scurve_ped[ch].Fill(da,dat); // fill the histogram for pedestal data
     386        }
     387        in_ped.close();
     388        in_ped.clear();
     389       
     390        //dac_ped[ch]  = hist_scurve_ped[ch].GetMaximumBin()+dac_min_ped; //absolute value
     391        dac_ped[ch]  = hist_scurve_ped[ch].GetBinCenter(hist_scurve_ped[ch].GetMaximumBin()); //absolute value
     392        dac_ped_sum += dac_ped[ch]; //absolute value
     393        //cerr << "dac_ped[" << ch << "] dac_ped_sum " << dac_ped[ch]
     394        //<< " " << dac_ped_sum << endl;
     395       
     396       
     397        hist_scurve[ch].Smooth(3); // Smooth the S curve
     398        Int_t binx=0;
     399        //hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
     400        //                                  dac_max-dac_min+1,kHalfMax/2.);
     401       
     402        // ??????
     403        hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
     404                                          dac_max-dac_min+1,2000);
     405        dac_half[ch] = hist_scurve[ch].GetBinCenter(binx); //absolute value
     406        if(dac_half[ch]<kXmin)dac_half[ch]=-1;
     407       
     408        // put the histograms of the smoothed S-Curves into two arrays (dac, data)
     409        //-------------------------------------------------------------------------
     410        Int_t i=0;
     411        for(i=0;i<kNdata;i++){
     412            dac[i]  = hist_scurve[ch].GetBinCenter(i+1);
     413            data[i] = hist_scurve[ch].GetBinContent(i+1);
     414        }
     415       
     416        /*** compute derive diff on S-Curve ***/
     417        Double_t diff[kNdata-1], dac_diff[kNdata-1];
     418        Int_t trigger_diff=0;
     419        for(i=0;i<kNdata-1;i++){
     420            diff[i]     = (data[i]-data[i+1])/step_dac;
     421            dac_diff[i] = dac[i]+step_dac/2;
     422            hist_diff[ch].Fill(dac_diff[i],diff[i]); // fills differentiated histograms of smoothed s-curves which is the single photoelectron spectrum
     423            if(dac_diff[i]>=fit_min){
     424                diff_sum += diff[i]*dac_diff[i];
     425                //cerr << "diff[" << i << "], dac_diff[" << i << "] "
     426                //<< diff[i] << " " << dac_diff[i] << endl;
     427                trigger_diff+=diff[i];
     428            }
     429        }
     430        //mean_diff[ch] = diff_sum/trigger_diff;
     431        hist_diff[ch].GetXaxis()->SetRangeUser(fit_min,fit_max);
     432        mean_diff[ch] = hist_diff[ch].GetMean();
     433        //cerr << "trigger_diff, mean_diff[" << ch << "] "
     434        //<< trigger_diff << " " << mean_diff[ch] << endl;
     435       
     436       
     437        /******************* GRAPHICS *****************/
     438       
     439       
     440        // Now start the analysis which consist in finding the peak in the differentiated spectrum
     441        c1->cd(ch+1);
     442        /*** fit scurves  ***/
     443        cerr << "ch " << ch;
     444       
     445        //    SDC 24/09/2013 :: Need to declare outside this variable in a compiled program
     446        TF1 *fit_pol ;
     447       
     448       
     449        // Fit a polynomial curve in a given range on the s-curve directly. This is a kind of DEBUG procedure.
     450       
     451        //if(ch==43 || ch==47){ // SDC: 25/09/2013 Why a special channel : the fit pol4 is performed for another DAC range
     452        /*   if(ch==170){
     453         fit_pol      = new TF1("fit_pol","pol4",150,fit_max_scurve);
     454         hist_scurve[ch].Fit(fit_pol,"","",150,fit_max_scurve);
     455         }else{ */
     456        fit_pol      = new TF1("fit_pol","pol4",fit_min_scurve,fit_max_scurve);
     457        hist_scurve[ch].Fit(fit_pol,"","",fit_min_scurve,fit_max_scurve);
     458        // }
     459       
     460       
     461       
     462        Double_t fit_par_pol4[5]; // container for fit parameter results
     463        fit_pol->GetParameters(fit_par_pol4);
     464        sfit_chiS[ch] =  fit_pol->GetChisquare(); // chi squared
     465        TF1 *func               = hist_scurve[ch].GetFunction("fit_pol");
     466        Double_t data_diff_fit;
     467        for(i=0;i<kNdataFit;i++){
     468            //if(ch==43 || ch==47){
     469            if(ch==170){
     470                dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+150-dac_min;
     471                data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+150-dac_min); // get the value of the fitted function
     472            }else{
     473                dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min;
     474                data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min); // get the value of the fitted function
     475                //cerr << "dac_fit[" << i << "], data_fit[" << i << "] "
     476                //<< dac_fit[i] << " " << data_fit[i] << endl;
     477            }
     478        }
     479       
     480        // derive the fitted function on the S-curve to obtain the SPE
     481        Double_t dac_diff_fit;
     482        for(i=0;i<kNdataFit-1;i++){
     483            data_diff_fit = (data_fit[i]-data_fit[i+1])/step_dac;
     484            hist_diff_fit[ch].Fill(dac_fit[i],data_diff_fit); // Fill differentiated fitted S-curve
     485        }
     486       
     487        TLine *lhalf_max   = new TLine(dac[0],kHalfMax,dac[kNdata-1],kHalfMax);
     488        TLine *ldac_half   = new TLine(dac_half[ch],0,dac_half[ch],
     489                                       (Double_t)kDrawMax);
     490        hist_scurve[ch].SetMaximum(kDrawMax);
     491        hist_scurve[ch].GetXaxis()->SetRangeUser(kDrawDacMin,kDrawDacMax);
     492        hist_scurve[ch].Draw();
     493        lhalf_max->Draw("same");
     494        ldac_half->Draw("same");
     495        if(dac_half[ch]>0){
     496            dac_half_sum += dac_half[ch]-dac_ped[ch]; //relative absolute value
     497            trigger++;
     498        }
     499       
     500        /*** Differential plot from scurve fitting function  ***/
     501        c5->cd(ch+1);
     502        hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve,fit_max_scurve-2*step_dac);
     503        Double_t dac_sfit_max=hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
     504        //if(dac_sfit_max<=fit_min_scurve+4*step_dac){
     505        if(dac_sfit_max<=fit_min_scurve+2*step_dac){
     506            hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve+kFitSadj,fit_max_scurve-2*step_dac);
     507            dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
     508            //if(dac_sfit[ch]<=fit_min_scurve+kFitSadj+4*step_dac)dac_sfit[ch]=-1;
     509            if(dac_sfit[ch]<fit_min_scurve+kFitSadj+step_dac){
     510                dac_sfit[ch]=-1;
     511                //dac_sfit[ch]=dac_half[ch];
     512                flag_fit[ch] = -1;
     513                cerr << "flag_fit[" << ch << "] " << flag_fit[ch]  << endl;
     514            }
     515        }else {
     516            dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
     517        }
     518       
     519        if(dac_sfit[ch]>fit_min_scurve+5){
     520            dac_sfit_sum += dac_sfit[ch]-dac_ped[ch];
     521            trigger_sfit++;
     522        }else{
     523            //dac_sfit[ch] = fit_min_scurve+10;
     524            dac_sfit[ch] = -1;
     525        }
     526        //cerr << "dac_sfit[" << ch << "] " << dac_sfit[ch] << endl;
     527        TLine *ldac_sfit      = new TLine(dac_sfit[ch],0,dac_sfit[ch],kDrawMaxFit);
     528        TLine *ldac_sfit_copy = new TLine(dac_sfit[ch],0,dac_sfit[ch],kDrawMax);
     529        ldac_sfit->SetLineColor(kMagenta);
     530        ldac_sfit_copy->SetLineColor(kMagenta);
     531        hist_diff_fit[ch].SetMaximum(kDrawMaxFit);
     532        hist_diff_fit[ch].SetMinimum(0);
     533        hist_diff_fit[ch].Draw();
     534        ldac_sfit->Draw("same");
     535       
     536        /**** differential plot ****/
     537        c2->cd(ch+1);
     538        //hist_diff[ch].Rebin(8);
     539        hist_diff[ch].Rebin(3);
     540        fit_sdiff     = new TF1("fit_sdiff","gaus",dac_min,dac_max);
     541        //fit_sdiff->SetLineWidth(1);
     542        //if(ch==8 || ch==18 || ch==44 || ch==56)
     543        mean_sum += mean_diff[ch]-dac_ped[ch];
     544        trigger_mean++;
     545       
     546        // fit the derivated polynom fit
     547        Double_t fit_par[3];
     548        //hist_diff[ch].Fit(fit_sdiff,"RQ0","",fit_min,fit_max);
     549        hist_diff[ch].Fit(fit_sdiff,"","",fit_min,fit_max);
     550        hist_diff[ch].GetXaxis()->SetRangeUser(kDrawDacMin,kDrawDacMax);
     551        hist_diff[ch].SetMinimum(0);
     552        hist_diff[ch].SetMaximum(kDrawMaxDiff);
     553        //hist_diff[ch].Rebin();
     554        hist_diff[ch].Draw();
     555        fit_sdiff->GetParameters(fit_par); // get the parameters for this fit
     556        dac_spe[ch] = fit_par[1];
     557        if(dac_spe[ch]>fit_min){
     558            dac_spe_sum += dac_spe[ch]-dac_ped[ch];
     559            //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
     560            trigger_fit++;
     561        }else{
     562            dac_spe[ch] = fit_min+10;
     563        }
     564        //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
     565        TLine *ldac_spe    = new TLine(dac_spe[ch],0,dac_spe[ch],
     566                                       (Double_t)kDrawMaxDiff);
     567        ldac_spe->SetLineColor(kRed);
     568        TLine *lmean_diff  = new TLine(mean_diff[ch],0,mean_diff[ch],
     569                                       (Double_t)kDrawMaxDiff);
     570        lmean_diff->SetLineColor(kCyan);
     571        TLine *lfit_min    = new TLine(fit_min,0,fit_min,kDrawMaxDiff);
     572        lfit_min->SetLineColor(kGray);
     573       
     574        ldac_spe->Draw();
     575        lmean_diff->Draw("same");
     576        lfit_min->Draw("same");
     577       
     578        c3->cd(1);
     579        if(ch==0){
     580            hist_scurve[ch].SetMaximum(kDrawMax);
     581            hist_scurve[ch].Draw();
     582        }else{
     583            /*if(ch==8 || ch==18 || ch==44 || ch==56)*/
     584            hist_scurve[ch].Draw("same");
     585        }
     586        if(ch==kDrawCh){
     587            c3->cd(2);
     588            hist_scurve[ch].SetMaximum(kDrawMax);
     589            hist_scurve[kDrawCh].Draw();
     590            lhalf_max->Draw("same");
     591            ldac_half->Draw("same");
     592            ldac_sfit_copy->Draw("same");
     593            c3->cd(3);
     594            hist_diff[kDrawCh].Draw();
     595            ldac_spe->Draw();
     596            lmean_diff->Draw("same");
     597            lfit_min->Draw("same");
     598            c3->cd(4);
     599            hist_diff_fit[kDrawCh].Draw();
     600            ldac_sfit->Draw("same");
     601        }
     602        /*
     603         if(ch<32){
     604         c6->cd(ch/4+1);
     605         hist_scurve[ch].SetMaximum(kDrawMax);
     606         if(ch%4==0){
     607         hist_scurve[ch].Draw();
     608         }else{
     609         hist_scurve[ch].Draw("same");
     610         }
     611         ldac_sfit_copy->Draw("same");
     612         }else{
     613         c7->cd(ch/4-7);
     614         hist_scurve[ch].SetMaximum(kDrawMax);
     615         if(ch%4==0){
     616         hist_scurve[ch].Draw();
     617         }else{
     618         hist_scurve[ch].Draw("same");
     619         }
     620         ldac_sfit_copy->Draw("same");
     621         }
     622         */
     623       
     624        c4->cd();
     625        if(ch==0){
     626            hist_scurve_ped[ch].SetMaximum(kDrawMaxPed);
     627            hist_scurve_ped[ch].Draw("l");
     628        }else{
     629            hist_scurve_ped[ch].Draw("lsame");
     630        }
     631       
     632    }// end of second loop on all the channels
     633   
     634    // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
     635   
     636    Float_t dac_ped_ave  ; //absolute value
     637    Float_t dac_half_ave  ;
     638    Float_t dac_spe_ave   ;
     639    Float_t dac_sfit_ave;
     640    Float_t mean_diff_ave = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
     641   
     642    if(!checkAve){
     643        dac_ped_ave   = dac_ped[kChRef]; //absolute value
     644        dac_half_ave  = dac_half[kChRef];
     645        dac_spe_ave   = dac_spe[kChRef];
     646        dac_sfit_ave  = dac_sfit[kChRef];
     647        //Float_t dac_sfit_ave  = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
     648        mean_diff_ave = mean_diff[kChRef];
     649    }else{
     650        dac_ped_ave  = dac_ped_sum/(Float_t)kNch; //absolute value
     651        dac_half_ave = dac_half_sum/(Float_t)trigger+dac_ped_ave;
     652        dac_spe_ave  = dac_spe_sum/(Float_t)trigger_fit+dac_ped_ave;
     653        mean_diff_ave = mean_sum/(Float_t)trigger_mean+dac_ped_ave;
     654        dac_sfit_ave = dac_sfit_sum/(Float_t)trigger_sfit+dac_ped_ave;
     655    }
     656   
     657   
     658    // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
     659   
     660    Double_t dac_half_spe;
     661   
     662    //trigger_mean=4;
     663    if(dac_half_spe_fix!=-1) {
     664        dac_half_spe = dac_half_spe_fix;
     665    } else {
     666        //Double_t dac_half_spe = dac_half_ave;
     667        dac_half_spe = dac_spe_ave;
     668    }
     669   
     670    TLine *ldac_ped_ave   = new TLine(dac_ped_ave,0,dac_ped_ave,
     671                                      (Double_t)kDrawMax);
     672    TLine *ldac_ped_ave_copy = new TLine(dac_ped_ave,0,dac_ped_ave,
     673                                         (Double_t)kDrawMaxDiff);
     674    TLine *ldac_half_ave  = new TLine(dac_half_ave,0,dac_half_ave,
     675                                      (Double_t)kDrawMax);
     676    TLine *ldac_half_spe  = new TLine(dac_half_spe,0,dac_half_spe,
     677                                      (Double_t)kDrawMax);
     678    TLine *lmean_diff_ave = new TLine(mean_diff_ave,0,mean_diff_ave,
     679                                      (Double_t)kDrawMax);
     680    TLine *ldac_sfit_ave  = new TLine(dac_sfit_ave,0,dac_sfit_ave,
     681                                      (Double_t)kDrawMax);
     682    ldac_ped_ave->SetLineColor(kOrange);
     683    ldac_half_ave->SetLineColor(kRed);
     684    lmean_diff_ave->SetLineColor(kCyan);
     685    //ldac_half_spe->SetLineColor(kGreen);
     686    ldac_sfit_ave->SetLineColor(kMagenta);
     687   
     688    c3->cd(1);
     689    ldac_ped_ave->Draw("same");
     690    ldac_half_ave->Draw("same");
     691    // lhalf_max->Draw("same");  SDC (24/09/2013) :: Not declared
     692    ldac_half_spe->Draw("same");
     693    lmean_diff_ave->Draw("same");
     694    ldac_sfit_ave->Draw("same");
     695   
     696    c3->cd(3);
     697    ldac_ped_ave_copy->SetLineColor(kOrange);
     698    ldac_ped_ave_copy->Draw("same");
     699   
     700    ofstream out(kFnameOut.Data());
     701   
     702    Int_t    gain_asic_diff, gain_asic_dacHalf, gain_asic_diff_mean, gain_asic_sfit;
     703    out << "title  ch  g_asic_dacHalf  g_asic_diff  g_asic_diff_mean  g_asic_sfit  dac_half[ch]  dac_ped[ch]  dac_spe[ch]  dac_half_ave dac_half_spe dac_ped_ave  dac_sfit_ave"
     704    << endl;
     705   
     706   
     707    // loop 3 on all channels to make final computations on the gain
     708    for(ch=0;ch<kNch;ch++){
     709        /*
     710         if(ch==8 || ch==44 || ch==56 || ch==18){
     711         gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
     712         (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
     713         }else{
     714         gain_asic_diff_mean=0;
     715         }
     716         */
     717        if(dac_half[ch]<0){
     718            gain_asic_dacHalf   = 0;
     719            gain_asic_diff      = 0;
     720            gain_asic_diff_mean = 0;
     721            gain_asic_sfit      = 0;
     722        }else{
     723            gain_asic_dacHalf  =(Int_t)((dac_half_ave-dac_ped_ave)/
     724                                        (dac_half[ch]-dac_ped[ch])*gain_org[ch]);
     725            gain_asic_diff     =(Int_t)((dac_half_spe-dac_ped_ave)/
     726                                        (dac_spe[ch]-dac_ped[ch])*gain_org[ch]);
     727            gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
     728                                        (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
     729            gain_asic_sfit     =(Int_t)((dac_sfit_ave-dac_ped_ave)/
     730                                        (dac_sfit[ch]-dac_ped[ch])*gain_org[ch]);
     731            //cerr << ch << " " << gain_org[ch] << endl;
     732        }
     733       
     734       
     735       
     736        out << "OUT_TABLE    "
     737        << ch << "    "
     738        << gain_asic_dacHalf   << "    "
     739        << gain_asic_diff      << "    "
     740        << gain_asic_diff_mean << "    "
     741        << gain_asic_sfit      << "    "
     742        << dac_half[ch]        << "    "
     743        << dac_ped[ch]         << "    "
     744        << (Int_t)dac_spe[ch]  << "    "
     745        << dac_half_ave        << "    "
     746        << dac_half_spe        << "    "
     747        << dac_ped_ave         << "    "
     748        << dac_sfit_ave        << "    "
     749        << sfit_chiS[ch]       << "    "
     750        << ch                  << "    "
     751        << flag_fit[ch]
     752        << endl;
     753        std::cout << "channel " << ch << " old gain " << gain_org[ch] << " new gain " << gain_asic_sfit << std::endl;
     754    }
     755   
     756   
     757   
     758    /*
     759     for(ch=kNch-1;ch>=0;ch--){
     760     cout << "OUT_TABLE_INV    "
     761     << (Int_t)(dac_half_ave/(dac_half[ch]-dac_ped)*64)
     762     << endl;
     763     }
     764     */
     765   
     766   
     767}
     768
     769//--------------------------------------------------------------------------------------------------------
     770
     771//--------------------------------------------------------------------------------------------------------
     772// Analysis Constructor
     773//--------------------------------------------------------------------------------------------------------
     774Analysis::Analysis():ch(0),dac_half_sum(0),dac_spe_sum(0),dac_sfit_sum(0),diff_sum(0),mean_sum(0),maxBinSum(0)
     775{
     776    if(DEBUG)
     777        cout << " >>>>>>> START Analysis Constructor " << endl;
     778   
     779    /*********************** SIGNAL RUN ********************/
     780    input             = new TString[kNch];
     781    input_ped         = new TString[kNch];
     782   
     783    // Init sums
     784    dac_half_sum      = 0;
     785    dac_ped_sum       = 0;
     786    dac_spe_sum       = 0;
     787    dac_sfit_sum      = 0;
     788    diff_sum          = 0;
     789    mean_sum          = 0;
     790    maxBinSum         = 0;
     791   
     792    if(DEBUG)
     793        cout << " >>>>>>> END Analysis Constructor " << endl;
     794}
     795//--------------------------------------------------------------------------------------------------------
     796
     797
     798//--------------------------------------------------------------------------------------------------------
     799// Analysis destructor
     800//--------------------------------------------------------------------------------------------------------
     801Analysis::~Analysis()
     802{
     803   
     804    if(DEBUG)
     805        cout << " >>>>>>> START Analysis Destructor " << endl;
     806   
     807    if(input!=0) delete input;
     808    if(input_ped!=0) delete input_ped;
     809   
     810    if(DEBUG)
     811        cout << " >>>>>>> END Analysis Destructor " << endl;
     812   
     813}
     814//--------------------------------------------------------------------------------------------------------
     815
     816
     817//--------------------------------------------------------------------------------------------------------
     818// Main public function of Analysis
     819//--------------------------------------------------------------------------------------------------------
     820/*
     821 void Analysis::PlotScurve(
     822 Int_t    kDrawCh=31,
     823 TString  fname=kFname,
     824 TString  fname_ped=kFnamePed, //ASIC_D ped_file
     825 Double_t dac_half_spe_fix = -1,
     826 TString  fname_gain="gain-org.txt",
     827 Int_t    fit_min=kFitMin,
     828 Int_t    fit_max=kFitMax,
     829 Int_t    fit_min_scurve=kFitMinS,
     830 Int_t    fit_max_scurve=kFitMaxS,
     831 Bool_t   kBoard=0, //0:ASIC_C, 1:ASIC_D
     832 Int_t    kChRef=kRefCh,//ASIC_C
     833 Bool_t   checkAve=0 //0:use ref pix, 1:check average
     834 )
     835 */
     836void Analysis:: PlotScurve(
     837                           Int_t    kDrawCh,
     838                           TString  fname,
     839                           TString  fname_ped, //ASIC_D ped_file
     840                           Double_t dac_half_spe_fix,
     841                           TString  fname_gain,
     842                           Int_t    fit_min,
     843                           Int_t    fit_max,
     844                           Int_t    fit_min_scurve,
     845                           Int_t    fit_max_scurve,
     846                           Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
     847                           Int_t    kChRef,//ASIC_C
     848                           Bool_t   checkAve //0:use ref pix, 1:check average
     849                           )
     850
     851
     852
     853{
     854   
     855    if(DEBUG)
     856    {
     857        cout << " >>>>>>> START Analysis::PlotScurve " << endl;
     858        cout << " \t - kDrawCh = " << kDrawCh << endl;
     859        cout << " \t - fname = " << fname << endl;
     860        cout << " \t - fname_ped = " << fname_ped <<endl;
     861        cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
     862       
     863        cout << " \t - fname_gain = " << fname_gain << endl;
     864        cout << " \t - fit_min = " << fit_min << endl;
     865        cout << " \t - fit_max = " << fit_max << endl;
     866        cout << " \t - fit_min_scurve = " << fit_min_scurve << endl;
     867        cout << " \t - fit_max_scurve = " << fit_max_scurve << endl;
     868        cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
     869        cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
     870        cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;
     871    }
     872    // reset all variables
     873    Init();
     874   
     875    //OpenGainFile(fname_gain);
     876   
     877    FillHistoName();
     878   
     879    ReadGainFile(fname_gain,kBoard);
     880   
     881    // Open the pedestal file and the S-Curve files for each channel and test if the file can be open
     882   
     883    trigger=0;
     884    trigger_fit=0;
     885    trigger_mean=0;
     886    trigger_sfit=0;
     887   
     888    // Big loop on all the channels : second loop
     889    //--------------------------------------------
     890    // build the input data file name for each channel fname-chXX.dat
     891    for(ch=0;ch<kNch;ch++)
     892    {
     893        ProcessSingleChannel(kDrawCh,fname,fname_ped,dac_half_spe_fix,fname_gain,
     894                             fit_min,fit_max,fit_min_scurve,fit_max_scurve,kBoard,kChRef);
     895    }
     896   
     897   
     898   
     899    ComputeDACAverage(dac_half_spe_fix ,kChRef,checkAve);
     900    OutputGain(kFnameOut);
     901   
     902    if(DEBUG)
     903        cout << " >>>>>>> END Analysis::PlotSCurve " << endl;
     904   
     905}
     906//--------------------------------------------------------------------------------------------------------
     907
     908
     909//--------------------------------------------------------------------------------------------------------
     910// Initialize the program by reseting counters or variables and create objects
     911//--------------------------------------------------------------------------------------------------------
     912void Analysis::Init()
     913{
     914   
     915    if(DEBUG)
     916        cout << " >>>>>>> START Analysis::Init " << endl;
     917   
     918    //if(input!=0) delete input;
     919    //if(input_ped!=0) delete input_ped;
     920   
     921    // recreate the name
     922    input             = new TString[kNch];
     923    input_ped         = new TString[kNch];
     924   
     925    // Init sums
     926    dac_half_sum      = 0;
     927    dac_ped_sum       = 0;
     928    dac_spe_sum       = 0;
     929    dac_sfit_sum      = 0;
     930    diff_sum          = 0;
     931    mean_sum          = 0;
     932    maxBinSum         = 0;
     933   
     934   
     935    // create canvas
     936   
     937    c1                = new TCanvas(kCanvasName1,kCanvasTitle1,200,20,1100,600); //smoothed S-curves
     938    c1->Divide(11,6);
     939    c2                = new TCanvas(kCanvasName2,kCanvasTitle2,300,50,1100,600); //differentiated histograms of smoothed s-curves
     940    c2->Divide(11,6);
     941    c5                = new TCanvas(kCanvasName5,kCanvasTitle5,400,70,1100,600);
     942    c5->Divide(11,6);
     943    c3                = new TCanvas(kCanvasName3,kCanvasTitle3,500,100,1200,300);
     944    c3->Divide(4,1);
     945    c4                = new TCanvas(kCanvasName4,kCanvasTitle4,900,400,300,300);
     946   
     947   
     948    // create histos and function
     949    hist_scurve       = new TH1D[kNch]; // smoothed S-curves
     950    hist_scurve_ped   = new TH1D[kNch]; // S-curves for pedestal runs
     951    hist_diff         = new TH1D[kNch]; // differentiated histograms of smoothed s-curves
     952    hist_diff_fit     = new TH1D[kNch]; // differentiated histograms of fitted s-curves
     953    fit_sdiff         = new TF1[kNch];  // function fitted on differentiated histograms of fitted s-curves
     954   
     955   
     956    for(ch=0;ch<kNch;ch++){
     957        dac_half[ch]=0;
     958        dac_ped[ch]=0;
     959        dac_spe[ch]=0;
     960        dac_sfit[ch]=0;
     961        gain_org[ch]=0;
     962        mean_diff[ch]=0;
     963        sfit_chiS[ch]=0;
     964        flag_fit[ch]=0;
     965    }
     966   
     967    if(DEBUG)
     968        cout << " >>>>>>> END Analysis::Init " << endl;
     969   
     970   
     971}
     972//--------------------------------------------------------------------------------------------------------
     973
     974//--------------------------------------------------------------------------------------------------------
     975// Initialize the program by reseting counters or variables and create objects
     976//--------------------------------------------------------------------------------------------------------
     977void Analysis::OpenGainFile(const TString fname_gain)
     978{
     979    if(DEBUG)
     980        cout << " >>>>>>> START Analysis::OpenGainFile " << endl;
     981   
     982   
     983   
     984    TString input_gain(fname_gain);
     985    ifstream in_gain(input_gain.Data());
     986    if(!in_gain.good()){
     987        cerr << "File " << input_gain << " open failed!" << endl;
     988        exit(-1);
     989    }
     990    in_gain.close();
     991    in_gain.clear();
     992   
     993    if(DEBUG)
     994        cout << " >>>>>>> END Analysis::OpenGainFile " << endl;
     995   
     996   
     997}
     998//--------------------------------------------------------------------------------------------------------
     999
     1000//--------------------------------------------------------------------------------------------------------
     1001// Fill Histogram name
     1002//--------------------------------------------------------------------------------------------------------
     1003void Analysis::FillHistoName()
     1004{
     1005   
     1006    if(DEBUG)
     1007        cout << " >>>>>>> START Analysis::FillHistoName " << endl;
     1008   
     1009   
     1010   
     1011    // first loop on channels to fill histonames
     1012    //----------------------
     1013    for(ch=0;ch<kNch;ch++){
     1014       
     1015        name_hist[ch]           ="hist";
     1016        name_hist[ch]          +=ch;
     1017        name_hist_ped[ch]       ="hist_ped";
     1018        name_hist_ped[ch]      +=ch;
     1019        name_hist_diff[ch]      ="hist_diff";
     1020        name_hist_diff[ch]     +=ch;
     1021        name_hist_diff_fit[ch]  ="hist_diff_fit";
     1022        name_hist_diff_fit[ch] +=ch;
     1023        name_fit_diff[ch]       ="fit_diff";
     1024        name_fit_diff[ch]      +=ch;
     1025    }
     1026   
     1027    if(DEBUG)
     1028        cout << " >>>>>>> END Analysis::FillHistoName " << endl;
     1029   
     1030   
     1031   
     1032}
     1033//--------------------------------------------------------------------------------------------------------
     1034
     1035
     1036//--------------------------------------------------------------------------------------------------------
     1037// Read the gain
     1038// input : the gain filename
     1039// output the gain_org[ch] array
     1040//--------------------------------------------------------------------------------------------------------
     1041
     1042void Analysis::ReadGainFile(const TString fname_gain,Bool_t kBoard)
     1043{
     1044   
     1045    if(DEBUG)
     1046    {
     1047        cout << " >>>>>>> START Analysis::ReadGainFile " << endl;
     1048        cout << " \t - fname_gain = " << fname_gain << endl;
     1049        cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
     1050       
     1051    }
     1052   
     1053   
     1054   
     1055    TString input_gain(fname_gain);
     1056    ifstream in_gain(input_gain.Data());
     1057    if(!in_gain.good()){
     1058        cerr << "File " << input_gain << " open failed!" << endl;
     1059        exit(-1);
     1060    }
     1061    // read gain
     1062    Double_t da;
     1063    Int_t count=0;
     1064    ch=0;
     1065    while(in_gain >> da){
     1066        if(!kBoard && (count>=128 && count<192)){
     1067            gain_org[ch]=da;
     1068            cerr << "!kBoard, count, gain_org[" << ch << "] "
     1069            << !kBoard << " " << count << " " << gain_org[ch] << endl;
     1070            ch++;
     1071        }
     1072        if(kBoard && (count>=192 && count<256)){
     1073            gain_org[ch]=da;
     1074            cerr << "kBoard, count, gain_org[" << ch << "] "
     1075            << kBoard << " " << count << " " << gain_org[ch] << endl;
     1076            ch++;
     1077        }
     1078        count++;
     1079    }// end if first loop on all channels
     1080    in_gain.close();
     1081    in_gain.clear();
     1082   
     1083    if(DEBUG)
     1084    {
     1085        cout << " >>>>>>> END Analysis::ReadGainFile " << endl;
     1086       
     1087    }
     1088   
     1089}
     1090//--------------------------------------------------------------------------------------------------------
     1091
     1092
     1093
     1094//--------------------------------------------------------------------------------------------------------
     1095// ProcessSingleChannel
     1096//--------------------------------------------------------------------------------------------------------
     1097void Analysis::ProcessSingleChannel(
     1098                                    Int_t    kDrawCh,
     1099                                    TString  fname,
     1100                                    TString  fname_ped, //ASIC_D ped_file
     1101                                    Double_t dac_half_spe_fix,
     1102                                    TString  fname_gain,
     1103                                    Int_t    fit_min,
     1104                                    Int_t    fit_max,
     1105                                    Int_t    fit_min_scurve,
     1106                                    Int_t    fit_max_scurve,
     1107                                    Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
     1108                                    Int_t    kChRef//ASIC_C
     1109                                    //  Bool_t   checkAve //0:use ref pix, 1:check average)
     1110                                    )
     1111{
     1112   
     1113   
     1114    if(DEBUG)
     1115    {
     1116        cout << " >>>>>>> START Analysis::ProcessSingleChannel ch=" << ch << endl;
     1117        cout << " \t - kDrawCh = " << kDrawCh << endl;
     1118        cout << " \t - fname = " << fname << endl;
     1119        cout << " \t - fname_ped = " << fname_ped <<endl;
     1120        cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
     1121       
     1122        cout << " \t - fname_gain = " << fname_gain << endl;
     1123        cout << " \t - fit_min = " << fit_min << endl;
     1124        cout << " \t - fit_max = " << fit_max << endl;
     1125        cout << " \t - fit_min_scurve = " << fit_min_scurve << endl;
     1126        cout << " \t - fit_max_scurve = " << fit_max_scurve << endl;
     1127        cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
     1128        cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
     1129        //    cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;L
     1130    }
     1131   
     1132   
     1133   
    2791134    input[ch] += fname;
    2801135    input[ch] += "-ch";
     
    2831138    ifstream in(input[ch].Data());
    2841139    if(!in.good()){
    285       cerr << "File " << input[ch] << " open failed!" << endl;
    286       return;
    287     }
    288 
     1140        cerr << "File " << input[ch] << " open failed!" << endl;
     1141        return;
     1142    }
     1143   
    2891144    // build the input data file name for each channel fname_ped-chXX.dat
    2901145    input_ped[ch] += fname_ped;
     
    2941149    ifstream in_ped(input_ped[ch].Data());
    2951150    if(!in_ped.good()){
    296       cerr << "File " << input_ped[ch] << " open failed!" << endl;
    297       return;
    298     }   
    299 
    300 
     1151        cerr << "File " << input_ped[ch] << " open failed!" << endl;
     1152        return;
     1153    }
     1154   
     1155   
    3011156   
    3021157    //10DAC~10mV
     
    3041159    // and count the number of events
    3051160    Double_t da, dat, dac_min, dac_max, step_dac;
    306     count = 0;
     1161    Int_t count = 0;
    3071162    while (in >> da >> dat ){
    308       if(count==0) dac_min  = da;
    309       if(count==1) step_dac = da - dac_min;
    310       dac_max=da;
    311       count++;
    312     }     
     1163        if(count==0) dac_min  = da;
     1164        if(count==1) step_dac = da - dac_min;
     1165        dac_max=da;
     1166        count++;
     1167    }
    3131168    in.close();
    3141169    in.clear();
    315 
    316     /*   
    317     if(ch==43 || ch==47){//for T1 ASIC_C
    318       cerr << "CHECK " << ch << " " << fit_min_scurve << endl;
    319       Int_t fit_min_scurve=140;
    320     } else {
    321       Int_t fit_min_scurve=kFitMinS;
    322     }
    323     */
    324 
     1170   
     1171    /*
     1172     if(ch==43 || ch==47){//for T1 ASIC_C
     1173     cerr << "CHECK " << ch << " " << fit_min_scurve << endl;
     1174     Int_t fit_min_scurve=140;
     1175     } else {
     1176     Int_t fit_min_scurve=kFitMinS;
     1177     }
     1178     */
     1179   
    3251180    // correct the scale for DAC range
    3261181    if(kXmin<dac_min)kXmin=dac_min;
     
    3291184    const Int_t kNdataFit = count;
    3301185    if(ch==0){
    331       cerr << "kNdata, dac_min, dac_max, step_dac: "
    332            << kNdata << " " << dac_min << " " << dac_max << " "
    333            << step_dac << endl;
     1186        cerr << "kNdata, dac_min, dac_max, step_dac: "
     1187        << kNdata << " " << dac_min << " " << dac_max << " "
     1188        << step_dac << endl;
    3341189    }
    3351190   
     
    3371192    count = 0;
    3381193    while (in_ped >> da >> dat ){
    339       if(count==0) dac_min_ped  = da;
    340       if(count==1) step_dac_ped = da - dac_min_ped;
    341       dac_max_ped=da;
    342       count++;
    343     }
    344      
     1194        if(count==0) dac_min_ped  = da;
     1195        if(count==1) step_dac_ped = da - dac_min_ped;
     1196        dac_max_ped=da;
     1197        count++;
     1198    }
     1199   
    3451200    in_ped.close();
    3461201    in_ped.clear();
    3471202   
    348 
     1203   
    3491204    const Int_t kNdataPed = count;
    3501205    if(ch==0){
    351       cerr << "kNdataPed, dac_min_ped, dac_max_ped, step_dac_ped: "
    352            << kNdataPed << " " << dac_min_ped << " " << dac_max_ped
    353            << " " << step_dac_ped << endl;
     1206        cerr << "kNdataPed, dac_min_ped, dac_max_ped, step_dac_ped: "
     1207        << kNdataPed << " " << dac_min_ped << " " << dac_max_ped
     1208        << " " << step_dac_ped << endl;
    3541209    }
    3551210   
     
    3571212    // SDC: 27/09/2013 hist_scurve is a pointer on an array of THD, thus the syntax is correct (no new operator)
    3581213   
    359 
     1214   
    3601215    hist_scurve[ch]      = TH1D(name_hist[ch],name_hist[ch],kNdata-1,
    361                                 dac_min,dac_max);
     1216                                dac_min,dac_max);
    3621217    hist_scurve_ped[ch]  = TH1D(name_hist_ped[ch],name_hist_ped[ch],
    363                                 kNdataPed-1,dac_min_ped,dac_max_ped);
     1218                                kNdataPed-1,dac_min_ped,dac_max_ped);
    3641219    hist_diff[ch]        = TH1D(name_hist_diff[ch],name_hist_diff[ch],
    365                                 kNdata-2,dac_min,dac_max);
     1220                                kNdata-2,dac_min,dac_max);
    3661221    hist_diff_fit[ch]    = TH1D(name_hist_diff_fit[ch],name_hist_diff_fit[ch],
    367                                 kNdataFit,fit_min_scurve,fit_max_scurve);
    368 
    369    
    370 
    371 
    372 
     1222                                kNdataFit,fit_min_scurve,fit_max_scurve);
     1223   
     1224   
     1225   
     1226   
     1227   
    3731228    // read the data for that channel for normal data
    3741229    in.open(input[ch].Data());
    3751230    Double_t       data[kNdata], dac[kNdata], dac_fit[kNdataFit], data_fit[kNdataFit];
    3761231    while (in >> da >> dat ){
    377       hist_scurve[ch].Fill(da,dat); // fill the histograms of the s-curve
     1232        hist_scurve[ch].Fill(da,dat); // fill the histograms of the s-curve
    3781233    }
    3791234    in.close();
    3801235    in.clear();
    381 
     1236   
    3821237    // read the data for pedestal data
    3831238    in_ped.open(input_ped[ch].Data());
    3841239    while (in_ped >> da >> dat ){
    385       hist_scurve_ped[ch].Fill(da,dat); // fill the histogram for pedestal data
     1240        hist_scurve_ped[ch].Fill(da,dat); // fill the histogram for pedestal data
    3861241    }
    3871242    in_ped.close();
    3881243    in_ped.clear();
    389 
     1244   
    3901245    //dac_ped[ch]  = hist_scurve_ped[ch].GetMaximumBin()+dac_min_ped; //absolute value
    3911246    dac_ped[ch]  = hist_scurve_ped[ch].GetBinCenter(hist_scurve_ped[ch].GetMaximumBin()); //absolute value
    3921247    dac_ped_sum += dac_ped[ch]; //absolute value
    393     //cerr << "dac_ped[" << ch << "] dac_ped_sum " << dac_ped[ch] 
     1248    //cerr << "dac_ped[" << ch << "] dac_ped_sum " << dac_ped[ch]
    3941249    //<< " " << dac_ped_sum << endl;
    3951250   
    396 
     1251   
    3971252    hist_scurve[ch].Smooth(3); // Smooth the S curve
    3981253    Int_t binx=0;
    3991254    //hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
    4001255    //                                  dac_max-dac_min+1,kHalfMax/2.);
    401 
     1256   
    4021257    // ??????
    4031258    hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
    404                                       dac_max-dac_min+1,2000);
     1259                                      dac_max-dac_min+1,2000);
    4051260    dac_half[ch] = hist_scurve[ch].GetBinCenter(binx); //absolute value
    4061261    if(dac_half[ch]<kXmin)dac_half[ch]=-1;
     
    4101265    Int_t i=0;
    4111266    for(i=0;i<kNdata;i++){
    412       dac[i]  = hist_scurve[ch].GetBinCenter(i+1);
    413       data[i] = hist_scurve[ch].GetBinContent(i+1);
    414     }
    415    
     1267        dac[i]  = hist_scurve[ch].GetBinCenter(i+1);
     1268        data[i] = hist_scurve[ch].GetBinContent(i+1);
     1269    }
     1270    
    4161271    /*** compute derive diff on S-Curve ***/
    4171272    Double_t diff[kNdata-1], dac_diff[kNdata-1];
    4181273    Int_t trigger_diff=0;
    4191274    for(i=0;i<kNdata-1;i++){
    420       diff[i]     = (data[i]-data[i+1])/step_dac;
    421       dac_diff[i] = dac[i]+step_dac/2;
    422       hist_diff[ch].Fill(dac_diff[i],diff[i]); // fills differentiated histograms of smoothed s-curves which is the single photoelectron spectrum
    423       if(dac_diff[i]>=fit_min){
    424         diff_sum += diff[i]*dac_diff[i];
    425         //cerr << "diff[" << i << "], dac_diff[" << i << "] "
    426         //<< diff[i] << " " << dac_diff[i] << endl;
    427         trigger_diff+=diff[i];
    428       }
     1275        diff[i]     = (data[i]-data[i+1])/step_dac;
     1276        dac_diff[i] = dac[i]+step_dac/2;
     1277        hist_diff[ch].Fill(dac_diff[i],diff[i]); // fills differentiated histograms of smoothed s-curves which is the single photoelectron spectrum
     1278        if(dac_diff[i]>=fit_min){
     1279            diff_sum += diff[i]*dac_diff[i];
     1280            //cerr << "diff[" << i << "], dac_diff[" << i << "] "
     1281            //<< diff[i] << " " << dac_diff[i] << endl;
     1282            trigger_diff+=diff[i];
     1283        }
    4291284    }
    4301285    //mean_diff[ch] = diff_sum/trigger_diff;
    4311286    hist_diff[ch].GetXaxis()->SetRangeUser(fit_min,fit_max);
    4321287    mean_diff[ch] = hist_diff[ch].GetMean();
    433     //cerr << "trigger_diff, mean_diff[" << ch << "] " 
     1288    //cerr << "trigger_diff, mean_diff[" << ch << "] "
    4341289    //<< trigger_diff << " " << mean_diff[ch] << endl;
    4351290   
    436    
     1291    
    4371292    /******************* GRAPHICS *****************/
    438 
    439 
     1293   
     1294   
    4401295    // Now start the analysis which consist in finding the peak in the differentiated spectrum
    441     c1->cd(ch+1);   
     1296    c1->cd(ch+1);
    4421297    /*** fit scurves  ***/
    4431298    cerr << "ch " << ch;
    444 
     1299   
    4451300    //    SDC 24/09/2013 :: Need to declare outside this variable in a compiled program
    4461301    TF1 *fit_pol ;
    447 
    448 
     1302   
     1303   
    4491304    // Fit a polynomial curve in a given range on the s-curve directly. This is a kind of DEBUG procedure.
    450 
     1305   
    4511306    //if(ch==43 || ch==47){ // SDC: 25/09/2013 Why a special channel : the fit pol4 is performed for another DAC range
    452  /*   if(ch==170){
    453       fit_pol      = new TF1("fit_pol","pol4",150,fit_max_scurve);
    454       hist_scurve[ch].Fit(fit_pol,"","",150,fit_max_scurve);
    455     }else{ */
    456       fit_pol      = new TF1("fit_pol","pol4",fit_min_scurve,fit_max_scurve);
    457       hist_scurve[ch].Fit(fit_pol,"","",fit_min_scurve,fit_max_scurve);
    458    // }
    459 
    460    
    461 
     1307    /*   if(ch==170){
     1308     fit_pol      = new TF1("fit_pol","pol4",150,fit_max_scurve);
     1309     hist_scurve[ch].Fit(fit_pol,"","",150,fit_max_scurve);
     1310     }else{ */
     1311    fit_pol      = new TF1("fit_pol","pol4",fit_min_scurve,fit_max_scurve);
     1312    hist_scurve[ch].Fit(fit_pol,"","",fit_min_scurve,fit_max_scurve);
     1313    // }
     1314   
     1315   
     1316   
    4621317    Double_t fit_par_pol4[5]; // container for fit parameter results
    4631318    fit_pol->GetParameters(fit_par_pol4);
     
    4661321    Double_t data_diff_fit;
    4671322    for(i=0;i<kNdataFit;i++){
    468       //if(ch==43 || ch==47){
    469       if(ch==170){
    470         dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+150-dac_min;
    471         data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+150-dac_min); // get the value of the fitted function
    472       }else{
    473         dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min;
    474         data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min); // get the value of the fitted function
    475         //cerr << "dac_fit[" << i << "], data_fit[" << i << "] "
    476         //<< dac_fit[i] << " " << data_fit[i] << endl;
    477       }
    478     }
    479 
     1323        //if(ch==43 || ch==47){
     1324        if(ch==170){
     1325            dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+150-dac_min;
     1326            data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+150-dac_min); // get the value of the fitted function
     1327        }else{
     1328            dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min;
     1329            data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min); // get the value of the fitted function
     1330            //cerr << "dac_fit[" << i << "], data_fit[" << i << "] "
     1331            //<< dac_fit[i] << " " << data_fit[i] << endl;
     1332        }
     1333    }
     1334   
    4801335    // derive the fitted function on the S-curve to obtain the SPE
    4811336    Double_t dac_diff_fit;
    4821337    for(i=0;i<kNdataFit-1;i++){
    483       data_diff_fit = (data_fit[i]-data_fit[i+1])/step_dac;
    484       hist_diff_fit[ch].Fill(dac_fit[i],data_diff_fit); // Fill differentiated fitted S-curve
    485     }
    486 
     1338        data_diff_fit = (data_fit[i]-data_fit[i+1])/step_dac;
     1339        hist_diff_fit[ch].Fill(dac_fit[i],data_diff_fit); // Fill differentiated fitted S-curve
     1340    }
     1341   
    4871342    TLine *lhalf_max   = new TLine(dac[0],kHalfMax,dac[kNdata-1],kHalfMax);
    4881343    TLine *ldac_half   = new TLine(dac_half[ch],0,dac_half[ch],
    489                                    (Double_t)kDrawMax);
     1344                                   (Double_t)kDrawMax);
    4901345    hist_scurve[ch].SetMaximum(kDrawMax);
    4911346    hist_scurve[ch].GetXaxis()->SetRangeUser(kDrawDacMin,kDrawDacMax);
    4921347    hist_scurve[ch].Draw();
    493     lhalf_max->Draw("same"); 
     1348    lhalf_max->Draw("same");
    4941349    ldac_half->Draw("same");
    4951350    if(dac_half[ch]>0){
    496       dac_half_sum += dac_half[ch]-dac_ped[ch]; //relative absolute value
    497       trigger++;
    498     }
    499 
     1351        dac_half_sum += dac_half[ch]-dac_ped[ch]; //relative absolute value
     1352        trigger++;
     1353    }
     1354   
    5001355    /*** Differential plot from scurve fitting function  ***/
    5011356    c5->cd(ch+1);
     
    5041359    //if(dac_sfit_max<=fit_min_scurve+4*step_dac){
    5051360    if(dac_sfit_max<=fit_min_scurve+2*step_dac){
    506       hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve+kFitSadj,fit_max_scurve-2*step_dac);
    507       dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
    508       //if(dac_sfit[ch]<=fit_min_scurve+kFitSadj+4*step_dac)dac_sfit[ch]=-1;
    509       if(dac_sfit[ch]<fit_min_scurve+kFitSadj+step_dac){
    510         dac_sfit[ch]=-1;
    511         //dac_sfit[ch]=dac_half[ch];
    512         flag_fit[ch] = -1;
    513         cerr << "flag_fit[" << ch << "] " << flag_fit[ch]  << endl;
    514       }
     1361        hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve+kFitSadj,fit_max_scurve-2*step_dac);
     1362        dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
     1363        //if(dac_sfit[ch]<=fit_min_scurve+kFitSadj+4*step_dac)dac_sfit[ch]=-1;
     1364        if(dac_sfit[ch]<fit_min_scurve+kFitSadj+step_dac){
     1365            dac_sfit[ch]=-1;
     1366            //dac_sfit[ch]=dac_half[ch];
     1367            flag_fit[ch] = -1;
     1368            cerr << "flag_fit[" << ch << "] " << flag_fit[ch]  << endl;
     1369        }
    5151370    }else {
    516       dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
    517     }
    518 
     1371        dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
     1372    }
     1373   
    5191374    if(dac_sfit[ch]>fit_min_scurve+5){
    520       dac_sfit_sum += dac_sfit[ch]-dac_ped[ch];
    521       trigger_sfit++;
     1375        dac_sfit_sum += dac_sfit[ch]-dac_ped[ch];
     1376        trigger_sfit++;
    5221377    }else{
    523       //dac_sfit[ch] = fit_min_scurve+10;
    524       dac_sfit[ch] = -1;
     1378        //dac_sfit[ch] = fit_min_scurve+10;
     1379        dac_sfit[ch] = -1;
    5251380    }
    5261381    //cerr << "dac_sfit[" << ch << "] " << dac_sfit[ch] << endl;
     
    5431398    mean_sum += mean_diff[ch]-dac_ped[ch];
    5441399    trigger_mean++;
    545 
    546     // fit the derivated polynom fit 
     1400   
     1401    // fit the derivated polynom fit
    5471402    Double_t fit_par[3];
    5481403    //hist_diff[ch].Fit(fit_sdiff,"RQ0","",fit_min,fit_max);
     
    5561411    dac_spe[ch] = fit_par[1];
    5571412    if(dac_spe[ch]>fit_min){
    558       dac_spe_sum += dac_spe[ch]-dac_ped[ch];
    559       //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
    560       trigger_fit++;
     1413        dac_spe_sum += dac_spe[ch]-dac_ped[ch];
     1414        //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
     1415        trigger_fit++;
    5611416    }else{
    562       dac_spe[ch] = fit_min+10;
     1417        dac_spe[ch] = fit_min+10;
    5631418    }
    5641419    //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
    5651420    TLine *ldac_spe    = new TLine(dac_spe[ch],0,dac_spe[ch],
    566                                 (Double_t)kDrawMaxDiff);
     1421                                  (Double_t)kDrawMaxDiff);
    5671422    ldac_spe->SetLineColor(kRed);
    5681423    TLine *lmean_diff  = new TLine(mean_diff[ch],0,mean_diff[ch],
    569                                    (Double_t)kDrawMaxDiff);
     1424                                   (Double_t)kDrawMaxDiff);
    5701425    lmean_diff->SetLineColor(kCyan);
    5711426    TLine *lfit_min    = new TLine(fit_min,0,fit_min,kDrawMaxDiff);
    5721427    lfit_min->SetLineColor(kGray);
    573 
     1428   
    5741429    ldac_spe->Draw();
    5751430    lmean_diff->Draw("same");
     
    5781433    c3->cd(1);
    5791434    if(ch==0){
    580       hist_scurve[ch].SetMaximum(kDrawMax);
    581       hist_scurve[ch].Draw();
     1435        hist_scurve[ch].SetMaximum(kDrawMax);
     1436        hist_scurve[ch].Draw();
    5821437    }else{
    583     /*if(ch==8 || ch==18 || ch==44 || ch==56)*/
    584       hist_scurve[ch].Draw("same");
     1438        /*if(ch==8 || ch==18 || ch==44 || ch==56)*/
     1439        hist_scurve[ch].Draw("same");
    5851440    }
    5861441    if(ch==kDrawCh){
    587       c3->cd(2);
    588       hist_scurve[ch].SetMaximum(kDrawMax);
    589       hist_scurve[kDrawCh].Draw();
    590       lhalf_max->Draw("same"); 
    591       ldac_half->Draw("same");
    592       ldac_sfit_copy->Draw("same");
    593       c3->cd(3);
    594       hist_diff[kDrawCh].Draw();     
    595       ldac_spe->Draw();
    596       lmean_diff->Draw("same");
    597       lfit_min->Draw("same");
    598       c3->cd(4);
    599       hist_diff_fit[kDrawCh].Draw();
    600        ldac_sfit->Draw("same");
    601       }
     1442        c3->cd(2);
     1443        hist_scurve[ch].SetMaximum(kDrawMax);
     1444        hist_scurve[kDrawCh].Draw();
     1445        lhalf_max->Draw("same");
     1446        ldac_half->Draw("same");
     1447        ldac_sfit_copy->Draw("same");
     1448        c3->cd(3);
     1449        hist_diff[kDrawCh].Draw();
     1450        ldac_spe->Draw();
     1451        lmean_diff->Draw("same");
     1452        lfit_min->Draw("same");
     1453        c3->cd(4);
     1454        hist_diff_fit[kDrawCh].Draw();
     1455        ldac_sfit->Draw("same");
     1456    }
    6021457    /*
    603     if(ch<32){
    604       c6->cd(ch/4+1);
    605       hist_scurve[ch].SetMaximum(kDrawMax);
    606       if(ch%4==0){
    607         hist_scurve[ch].Draw();
    608       }else{
    609         hist_scurve[ch].Draw("same");
    610       }
    611       ldac_sfit_copy->Draw("same");
    612     }else{
    613       c7->cd(ch/4-7);
    614       hist_scurve[ch].SetMaximum(kDrawMax);
    615       if(ch%4==0){
    616         hist_scurve[ch].Draw();
    617       }else{
    618         hist_scurve[ch].Draw("same");
    619       }
    620       ldac_sfit_copy->Draw("same");
    621       }
    622     */
     1458     if(ch<32){
     1459     c6->cd(ch/4+1);
     1460     hist_scurve[ch].SetMaximum(kDrawMax);
     1461     if(ch%4==0){
     1462     hist_scurve[ch].Draw();
     1463     }else{
     1464     hist_scurve[ch].Draw("same");
     1465     }
     1466     ldac_sfit_copy->Draw("same");
     1467     }else{
     1468     c7->cd(ch/4-7);
     1469     hist_scurve[ch].SetMaximum(kDrawMax);
     1470     if(ch%4==0){
     1471     hist_scurve[ch].Draw();
     1472     }else{
     1473     hist_scurve[ch].Draw("same");
     1474     }
     1475     ldac_sfit_copy->Draw("same");
     1476     }
     1477     */
    6231478   
    6241479    c4->cd();
    6251480    if(ch==0){
    626       hist_scurve_ped[ch].SetMaximum(kDrawMaxPed);
    627       hist_scurve_ped[ch].Draw("l");
     1481        hist_scurve_ped[ch].SetMaximum(kDrawMaxPed);
     1482        hist_scurve_ped[ch].Draw("l");
    6281483    }else{
    629       hist_scurve_ped[ch].Draw("lsame");
    630     }
    631    
    632   }// end of second loop on all the channels
    633 
    634   // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
    635 
    636   Float_t dac_ped_ave  ; //absolute value
    637   Float_t dac_half_ave  ;
    638   Float_t dac_spe_ave   ;
    639   Float_t dac_sfit_ave;
    640   Float_t mean_diff_ave = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
    641  
    642   if(!checkAve){
    643     dac_ped_ave   = dac_ped[kChRef]; //absolute value
    644     dac_half_ave  = dac_half[kChRef];
    645     dac_spe_ave   = dac_spe[kChRef];
    646     dac_sfit_ave  = dac_sfit[kChRef];
    647     //Float_t dac_sfit_ave  = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
    648     mean_diff_ave = mean_diff[kChRef];
    649   }else{
    650     dac_ped_ave  = dac_ped_sum/(Float_t)kNch; //absolute value
    651     dac_half_ave = dac_half_sum/(Float_t)trigger+dac_ped_ave;
    652     dac_spe_ave  = dac_spe_sum/(Float_t)trigger_fit+dac_ped_ave;
    653     mean_diff_ave = mean_sum/(Float_t)trigger_mean+dac_ped_ave;
    654     dac_sfit_ave = dac_sfit_sum/(Float_t)trigger_sfit+dac_ped_ave;
    655   }
    656 
    657 
    658 // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
    659 
    660   Double_t dac_half_spe;
    661 
    662   //trigger_mean=4;
    663   if(dac_half_spe_fix!=-1) {
    664     dac_half_spe = dac_half_spe_fix;
    665   } else {
    666     //Double_t dac_half_spe = dac_half_ave;
    667     dac_half_spe = dac_spe_ave;
    668   }
    669 
    670   TLine *ldac_ped_ave   = new TLine(dac_ped_ave,0,dac_ped_ave,
    671                                    (Double_t)kDrawMax);
    672   TLine *ldac_ped_ave_copy = new TLine(dac_ped_ave,0,dac_ped_ave,
    673                                    (Double_t)kDrawMaxDiff);
    674   TLine *ldac_half_ave  = new TLine(dac_half_ave,0,dac_half_ave,
    675                                    (Double_t)kDrawMax);
    676   TLine *ldac_half_spe  = new TLine(dac_half_spe,0,dac_half_spe,
    677                                    (Double_t)kDrawMax);
    678   TLine *lmean_diff_ave = new TLine(mean_diff_ave,0,mean_diff_ave,
    679                                    (Double_t)kDrawMax);
    680   TLine *ldac_sfit_ave  = new TLine(dac_sfit_ave,0,dac_sfit_ave,
    681                                    (Double_t)kDrawMax);
    682   ldac_ped_ave->SetLineColor(kOrange);
    683   ldac_half_ave->SetLineColor(kRed);
    684   lmean_diff_ave->SetLineColor(kCyan);
    685   //ldac_half_spe->SetLineColor(kGreen);
    686   ldac_sfit_ave->SetLineColor(kMagenta); 
    687 
    688   c3->cd(1);
    689   ldac_ped_ave->Draw("same");
    690   ldac_half_ave->Draw("same");
    691   // lhalf_max->Draw("same");  SDC (24/09/2013) :: Not declared
    692   ldac_half_spe->Draw("same");
    693   lmean_diff_ave->Draw("same");
    694   ldac_sfit_ave->Draw("same");
    695 
    696   c3->cd(3);
    697   ldac_ped_ave_copy->SetLineColor(kOrange);
    698   ldac_ped_ave_copy->Draw("same");
    699 
    700   ofstream out(kFnameOut.Data());
    701 
    702   Int_t    gain_asic_diff, gain_asic_dacHalf, gain_asic_diff_mean, gain_asic_sfit;
    703   out << "title  ch  g_asic_dacHalf  g_asic_diff  g_asic_diff_mean  g_asic_sfit  dac_half[ch]  dac_ped[ch]  dac_spe[ch]  dac_half_ave dac_half_spe dac_ped_ave  dac_sfit_ave"
    704        << endl;
    705 
    706 
    707   // loop 3 on all channels to make final computations on the gain
    708   for(ch=0;ch<kNch;ch++){
    709     /*
    710     if(ch==8 || ch==44 || ch==56 || ch==18){
    711       gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
    712                                   (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
    713     }else{
    714       gain_asic_diff_mean=0;
    715     }
    716     */
    717     if(dac_half[ch]<0){
    718       gain_asic_dacHalf   = 0;
    719       gain_asic_diff      = 0;
    720       gain_asic_diff_mean = 0;
    721       gain_asic_sfit      = 0;
    722     }else{
    723       gain_asic_dacHalf  =(Int_t)((dac_half_ave-dac_ped_ave)/
    724                                 (dac_half[ch]-dac_ped[ch])*gain_org[ch]);
    725       gain_asic_diff     =(Int_t)((dac_half_spe-dac_ped_ave)/
    726                              (dac_spe[ch]-dac_ped[ch])*gain_org[ch]);
    727       gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
    728                                   (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
    729       gain_asic_sfit     =(Int_t)((dac_sfit_ave-dac_ped_ave)/
    730                              (dac_sfit[ch]-dac_ped[ch])*gain_org[ch]);
    731       //cerr << ch << " " << gain_org[ch] << endl;
    732     }
    733 
    734    
    735 
    736     out << "OUT_TABLE    "
    737          << ch << "    "
    738          << gain_asic_dacHalf   << "    "
    739          << gain_asic_diff      << "    "
    740          << gain_asic_diff_mean << "    "
    741          << gain_asic_sfit      << "    "
    742          << dac_half[ch]        << "    "
    743          << dac_ped[ch]         << "    "
    744          << (Int_t)dac_spe[ch]  << "    "
    745          << dac_half_ave        << "    "
    746          << dac_half_spe        << "    "
    747          << dac_ped_ave         << "    "
    748          << dac_sfit_ave        << "    "
    749          << sfit_chiS[ch]       << "    "
    750          << ch                  << "    "
    751          << flag_fit[ch]       
    752          << endl;
    753       std::cout << "channel " << ch << " old gain " << gain_org[ch] << " new gain " << gain_asic_sfit << std::endl;
    754   }
    755 
    756 
    757 
    758   /*
    759     for(ch=kNch-1;ch>=0;ch--){
    760     cout << "OUT_TABLE_INV    "
    761     << (Int_t)(dac_half_ave/(dac_half[ch]-dac_ped)*64)
    762       << endl;
    763       }
    764   */
    765  
    766  
    767 }
    768 
    769 //--------------------------------------------------------------------------------------------------------
    770 
    771 //--------------------------------------------------------------------------------------------------------
    772 // Analysis Constructor
    773 //--------------------------------------------------------------------------------------------------------
    774 Analysis::Analysis():ch(0),dac_half_sum(0),dac_spe_sum(0),dac_sfit_sum(0),diff_sum(0),mean_sum(0),maxBinSum(0)
    775 {
    776   if(DEBUG)
    777     cout << " >>>>>>> START Analysis Constructor " << endl;
    778 
    779  /*********************** SIGNAL RUN ********************/
    780   input             = new TString[kNch];
    781   input_ped         = new TString[kNch];
    782 
    783   // Init sums
    784   dac_half_sum      = 0;
    785   dac_ped_sum       = 0;
    786   dac_spe_sum       = 0;
    787   dac_sfit_sum      = 0;
    788   diff_sum          = 0;
    789   mean_sum          = 0;
    790   maxBinSum         = 0;
    791 
    792  if(DEBUG)
    793     cout << " >>>>>>> END Analysis Constructor " << endl;
    794 }
    795 //--------------------------------------------------------------------------------------------------------
    796 
    797 
    798 //--------------------------------------------------------------------------------------------------------
    799 // Analysis destructor
    800 //--------------------------------------------------------------------------------------------------------
    801 Analysis::~Analysis()
    802 {
    803  
    804  if(DEBUG)
    805     cout << " >>>>>>> START Analysis Destructor " << endl;
    806 
    807   if(input!=0) delete input;
    808   if(input_ped!=0) delete input_ped;
    809 
    810 if(DEBUG)
    811     cout << " >>>>>>> END Analysis Destructor " << endl;
    812 
    813 }
    814 //--------------------------------------------------------------------------------------------------------
    815 
    816 
    817 //--------------------------------------------------------------------------------------------------------
    818 // Main public function of Analysis
    819 //--------------------------------------------------------------------------------------------------------
    820 /*
    821 void Analysis::PlotScurve(
    822                 Int_t    kDrawCh=31,
    823                 TString  fname=kFname,
    824                 TString  fname_ped=kFnamePed, //ASIC_D ped_file
    825                 Double_t dac_half_spe_fix = -1,
    826                 TString  fname_gain="gain-org.txt",
    827                 Int_t    fit_min=kFitMin,
    828                 Int_t    fit_max=kFitMax,
    829                 Int_t    fit_min_scurve=kFitMinS,
    830                 Int_t    fit_max_scurve=kFitMaxS,
    831                 Bool_t   kBoard=0, //0:ASIC_C, 1:ASIC_D
    832                 Int_t    kChRef=kRefCh,//ASIC_C
    833                 Bool_t   checkAve=0 //0:use ref pix, 1:check average
    834                 )
    835 */
    836 void Analysis:: PlotScurve(
    837                 Int_t    kDrawCh,
    838                 TString  fname,
    839                 TString  fname_ped, //ASIC_D ped_file
    840                 Double_t dac_half_spe_fix,
    841                 TString  fname_gain,
    842                 Int_t    fit_min,
    843                 Int_t    fit_max,
    844                 Int_t    fit_min_scurve,
    845                 Int_t    fit_max_scurve,
    846                 Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
    847                 Int_t    kChRef,//ASIC_C
    848                 Bool_t   checkAve //0:use ref pix, 1:check average
    849                 )
    850 
    851 
    852 
    853 {
    854 
    855 if(DEBUG)
    856 {
    857     cout << " >>>>>>> START Analysis::PlotScurve " << endl;
    858     cout << " \t - kDrawCh = " << kDrawCh << endl;
    859     cout << " \t - fname = " << fname << endl;
    860     cout << " \t - fname_ped = " << fname_ped <<endl;
    861     cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
    862 
    863     cout << " \t - fname_gain = " << fname_gain << endl;
    864     cout << " \t - fit_min = " << fit_min << endl;
    865     cout << " \t - fit_max = " << fit_max << endl;
    866     cout << " \t - fit_min_scurve = " << fit_min_scurve << endl;
    867     cout << " \t - fit_max_scurve = " << fit_max_scurve << endl;
    868     cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
    869     cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
    870     cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;     
    871 }
    872   // reset all variables
    873   Init();
    874 
    875   //OpenGainFile(fname_gain);
    876 
    877   FillHistoName();
    878 
    879   ReadGainFile(fname_gain,kBoard);
    880 
    881 // Open the pedestal file and the S-Curve files for each channel and test if the file can be open
    882 
    883   trigger=0;
    884   trigger_fit=0;
    885   trigger_mean=0;
    886   trigger_sfit=0;
    887 
    888   // Big loop on all the channels : second loop
    889   //--------------------------------------------
    890   // build the input data file name for each channel fname-chXX.dat
    891   for(ch=0;ch<kNch;ch++)
     1484        hist_scurve_ped[ch].Draw("lsame");
     1485    }
     1486   
     1487   
     1488    if(DEBUG)
    8921489    {
    893       ProcessSingleChannel(kDrawCh,fname,fname_ped,dac_half_spe_fix,fname_gain,
    894                            fit_min,fit_max,fit_min_scurve,fit_max_scurve,kBoard,kChRef);
    895       }
    896 
    897 
    898  
    899   ComputeDACAverage(dac_half_spe_fix ,kChRef,checkAve);
    900   OutputGain(kFnameOut);
    901 
    902 if(DEBUG)
    903     cout << " >>>>>>> END Analysis::PlotSCurve " << endl;
    904 
    905 }
    906 //--------------------------------------------------------------------------------------------------------
    907 
    908 
    909 //--------------------------------------------------------------------------------------------------------
    910 // Initialize the program by reseting counters or variables and create objects
    911 //--------------------------------------------------------------------------------------------------------
    912 void Analysis::Init()
    913 {
    914 
    915 if(DEBUG)
    916     cout << " >>>>>>> START Analysis::Init " << endl;
    917 
    918 //if(input!=0) delete input;
    919 //if(input_ped!=0) delete input_ped;
    920 
    921   // recreate the name
    922   input             = new TString[kNch];
    923   input_ped         = new TString[kNch];
    924 
    925   // Init sums
    926   dac_half_sum      = 0;
    927   dac_ped_sum       = 0;
    928   dac_spe_sum       = 0;
    929   dac_sfit_sum      = 0;
    930   diff_sum          = 0;
    931   mean_sum          = 0;
    932   maxBinSum         = 0;
    933 
    934 
    935   // create canvas
    936 
    937   c1                = new TCanvas(kCanvasName1,kCanvasTitle1,200,20,1100,600); //smoothed S-curves
    938   c1->Divide(11,6);
    939   c2                = new TCanvas(kCanvasName2,kCanvasTitle2,300,50,1100,600); //differentiated histograms of smoothed s-curves
    940   c2->Divide(11,6);
    941   c5                = new TCanvas(kCanvasName5,kCanvasTitle5,400,70,1100,600);
    942   c5->Divide(11,6);
    943   c3                = new TCanvas(kCanvasName3,kCanvasTitle3,500,100,1200,300);
    944   c3->Divide(4,1);
    945   c4                = new TCanvas(kCanvasName4,kCanvasTitle4,900,400,300,300);
    946 
    947 
    948   // create histos and function
    949   hist_scurve       = new TH1D[kNch]; // smoothed S-curves
    950   hist_scurve_ped   = new TH1D[kNch]; // S-curves for pedestal runs
    951   hist_diff         = new TH1D[kNch]; // differentiated histograms of smoothed s-curves
    952   hist_diff_fit     = new TH1D[kNch]; // differentiated histograms of fitted s-curves
    953   fit_sdiff         = new TF1[kNch];  // function fitted on differentiated histograms of fitted s-curves
    954 
    955 
    956  for(ch=0;ch<kNch;ch++){
    957     dac_half[ch]=0;
    958     dac_ped[ch]=0;
    959     dac_spe[ch]=0;
    960     dac_sfit[ch]=0;
    961     gain_org[ch]=0;
    962     mean_diff[ch]=0;
    963     sfit_chiS[ch]=0;
    964     flag_fit[ch]=0;
    965  }
    966 
    967 if(DEBUG)
    968     cout << " >>>>>>> END Analysis::Init " << endl;
    969 
    970  
    971 }
    972 //--------------------------------------------------------------------------------------------------------
    973 
    974 //--------------------------------------------------------------------------------------------------------
    975 // Initialize the program by reseting counters or variables and create objects
    976 //--------------------------------------------------------------------------------------------------------
    977 void Analysis::OpenGainFile(const TString fname_gain)
    978 {
    979 if(DEBUG)
    980     cout << " >>>>>>> START Analysis::OpenGainFile " << endl;
    981 
    982 
    983 
    984   TString input_gain(fname_gain);
    985   ifstream in_gain(input_gain.Data());
    986   if(!in_gain.good()){
    987     cerr << "File " << input_gain << " open failed!" << endl;
    988     exit(-1);
    989   }
    990   in_gain.close();
    991   in_gain.clear();
    992 
    993 if(DEBUG)
    994     cout << " >>>>>>> END Analysis::OpenGainFile " << endl;
    995 
    996 
    997 }
    998 //--------------------------------------------------------------------------------------------------------
    999 
    1000 //--------------------------------------------------------------------------------------------------------
    1001 // Fill Histogram name
    1002 //--------------------------------------------------------------------------------------------------------
    1003 void Analysis::FillHistoName()
    1004 {
    1005 
    1006 if(DEBUG)
    1007     cout << " >>>>>>> START Analysis::FillHistoName " << endl;
    1008 
    1009 
    1010  
    1011 // first loop on channels to fill histonames
    1012   //----------------------
    1013   for(ch=0;ch<kNch;ch++){
    1014    
    1015     name_hist[ch]           ="hist";
    1016     name_hist[ch]          +=ch;
    1017     name_hist_ped[ch]       ="hist_ped";
    1018     name_hist_ped[ch]      +=ch;
    1019     name_hist_diff[ch]      ="hist_diff";
    1020     name_hist_diff[ch]     +=ch;
    1021     name_hist_diff_fit[ch]  ="hist_diff_fit";
    1022     name_hist_diff_fit[ch] +=ch;
    1023     name_fit_diff[ch]       ="fit_diff";
    1024     name_fit_diff[ch]      +=ch;
    1025   }
    1026 
    1027 if(DEBUG)
    1028     cout << " >>>>>>> END Analysis::FillHistoName " << endl;
    1029 
    1030 
    1031 
    1032 }
    1033 //--------------------------------------------------------------------------------------------------------
    1034 
    1035 
    1036 //--------------------------------------------------------------------------------------------------------
    1037 // Read the gain
    1038 // input : the gain filename
    1039 // output the gain_org[ch] array
    1040 //--------------------------------------------------------------------------------------------------------
    1041 
    1042 void Analysis::ReadGainFile(const TString fname_gain,Bool_t kBoard)
    1043 {
    1044 
    1045 if(DEBUG)
    1046 {
    1047     cout << " >>>>>>> START Analysis::ReadGainFile " << endl;
    1048     cout << " \t - fname_gain = " << fname_gain << endl;
    1049     cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
    1050        
    1051 }
    1052 
    1053 
    1054 
    1055   TString input_gain(fname_gain);
    1056   ifstream in_gain(input_gain.Data());
    1057   if(!in_gain.good()){
    1058     cerr << "File " << input_gain << " open failed!" << endl;
    1059     exit(-1);
    1060   }
    1061  // read gain
    1062   Double_t da;
    1063   Int_t count=0;
    1064   ch=0;
    1065   while(in_gain >> da){
    1066     if(!kBoard && (count>=128 && count<192)){
    1067       gain_org[ch]=da;
    1068       cerr << "!kBoard, count, gain_org[" << ch << "] "
    1069            << !kBoard << " " << count << " " << gain_org[ch] << endl;
    1070       ch++;
    1071     }
    1072     if(kBoard && (count>=192 && count<256)){
    1073       gain_org[ch]=da;
    1074       cerr << "kBoard, count, gain_org[" << ch << "] "
    1075            << kBoard << " " << count << " " << gain_org[ch] << endl;
    1076       ch++;
    1077     }
    1078     count++;
    1079   }// end if first loop on all channels     
    1080   in_gain.close();
    1081   in_gain.clear();
    1082  
    1083 if(DEBUG)
    1084 {
    1085     cout << " >>>>>>> END Analysis::ReadGainFile " << endl;
    1086        
    1087 }
    1088 
    1089 }
    1090 //--------------------------------------------------------------------------------------------------------
    1091 
    1092 
    1093 
    1094 //--------------------------------------------------------------------------------------------------------
    1095 // ProcessSingleChannel
    1096 //--------------------------------------------------------------------------------------------------------
    1097 void Analysis::ProcessSingleChannel(
    1098                 Int_t    kDrawCh,
    1099                 TString  fname,
    1100                 TString  fname_ped, //ASIC_D ped_file
    1101                 Double_t dac_half_spe_fix,
    1102                 TString  fname_gain,
    1103                 Int_t    fit_min,
    1104                 Int_t    fit_max,
    1105                 Int_t    fit_min_scurve,
    1106                 Int_t    fit_max_scurve,
    1107                 Bool_t   kBoard, //0:ASIC_C, 1:ASIC_D
    1108                 Int_t    kChRef//ASIC_C
    1109                 //      Bool_t   checkAve //0:use ref pix, 1:check average)
    1110                                     )
    1111 {
    1112 
    1113 
    1114 if(DEBUG)
    1115 {
    1116   cout << " >>>>>>> START Analysis::ProcessSingleChannel ch=" << ch << endl;
    1117     cout << " \t - kDrawCh = " << kDrawCh << endl;
    1118     cout << " \t - fname = " << fname << endl;
    1119     cout << " \t - fname_ped = " << fname_ped <<endl;
    1120     cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
    1121 
    1122     cout << " \t - fname_gain = " << fname_gain << endl;
    1123     cout << " \t - fit_min = " << fit_min << endl;
    1124     cout << " \t - fit_max = " << fit_max << endl;
    1125     cout << " \t - fit_min_scurve = " << fit_min_scurve << endl;
    1126     cout << " \t - fit_max_scurve = " << fit_max_scurve << endl;
    1127     cout << " \t - kBoard = " << kBoard << " 0:ASIC_C, 1:ASIC_D " << endl;
    1128     cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
    1129     //    cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;L       
    1130 }
    1131 
    1132 
    1133 
    1134     input[ch] += fname;
    1135     input[ch] += "-ch";
    1136     input[ch] += ch;
    1137     input[ch] += ".dat";
    1138     ifstream in(input[ch].Data());
    1139     if(!in.good()){
    1140       cerr << "File " << input[ch] << " open failed!" << endl;
    1141       return;
    1142     }
    1143 
    1144     // build the input data file name for each channel fname_ped-chXX.dat
    1145     input_ped[ch] += fname_ped;
    1146     input_ped[ch] += "-ch";
    1147     input_ped[ch] += ch;
    1148     input_ped[ch] += ".dat";
    1149     ifstream in_ped(input_ped[ch].Data());
    1150     if(!in_ped.good()){
    1151       cerr << "File " << input_ped[ch] << " open failed!" << endl;
    1152       return;
    1153     }   
    1154 
    1155 
    1156    
    1157     //10DAC~10mV
    1158     // determine the DAC range in dac_min, dac_max, step_dac from file input_ped
    1159     // and count the number of events
    1160     Double_t da, dat, dac_min, dac_max, step_dac;
    1161     Int_t count = 0;
    1162     while (in >> da >> dat ){
    1163       if(count==0) dac_min  = da;
    1164       if(count==1) step_dac = da - dac_min;
    1165       dac_max=da;
    1166       count++;
    1167     }     
    1168     in.close();
    1169     in.clear();
    1170 
    1171     /*   
    1172     if(ch==43 || ch==47){//for T1 ASIC_C
    1173       cerr << "CHECK " << ch << " " << fit_min_scurve << endl;
    1174       Int_t fit_min_scurve=140;
    1175     } else {
    1176       Int_t fit_min_scurve=kFitMinS;
    1177     }
    1178     */
    1179 
    1180     // correct the scale for DAC range
    1181     if(kXmin<dac_min)kXmin=dac_min;
    1182     const Int_t kNdata    = count;
    1183     count = (fit_max_scurve-fit_min_scurve)/step_dac;
    1184     const Int_t kNdataFit = count;
    1185     if(ch==0){
    1186       cerr << "kNdata, dac_min, dac_max, step_dac: "
    1187            << kNdata << " " << dac_min << " " << dac_max << " "
    1188            << step_dac << endl;
    1189     }
    1190    
    1191     Double_t dac_min_ped, dac_max_ped, step_dac_ped;
    1192     count = 0;
    1193     while (in_ped >> da >> dat ){
    1194       if(count==0) dac_min_ped  = da;
    1195       if(count==1) step_dac_ped = da - dac_min_ped;
    1196       dac_max_ped=da;
    1197       count++;
    1198     }
    1199      
    1200     in_ped.close();
    1201     in_ped.clear();
    1202    
    1203 
    1204     const Int_t kNdataPed = count;
    1205     if(ch==0){
    1206       cerr << "kNdataPed, dac_min_ped, dac_max_ped, step_dac_ped: "
    1207            << kNdataPed << " " << dac_min_ped << " " << dac_max_ped
    1208            << " " << step_dac_ped << endl;
    1209     }
    1210    
    1211     // book histograms
    1212     // SDC: 27/09/2013 hist_scurve is a pointer on an array of THD, thus the syntax is correct (no new operator)
    1213    
    1214 
    1215     hist_scurve[ch]      = TH1D(name_hist[ch],name_hist[ch],kNdata-1,
    1216                                 dac_min,dac_max);
    1217     hist_scurve_ped[ch]  = TH1D(name_hist_ped[ch],name_hist_ped[ch],
    1218                                 kNdataPed-1,dac_min_ped,dac_max_ped);
    1219     hist_diff[ch]        = TH1D(name_hist_diff[ch],name_hist_diff[ch],
    1220                                 kNdata-2,dac_min,dac_max);
    1221     hist_diff_fit[ch]    = TH1D(name_hist_diff_fit[ch],name_hist_diff_fit[ch],
    1222                                 kNdataFit,fit_min_scurve,fit_max_scurve);
    1223 
    1224    
    1225 
    1226 
    1227 
    1228     // read the data for that channel for normal data
    1229     in.open(input[ch].Data());
    1230     Double_t       data[kNdata], dac[kNdata], dac_fit[kNdataFit], data_fit[kNdataFit];
    1231     while (in >> da >> dat ){
    1232       hist_scurve[ch].Fill(da,dat); // fill the histograms of the s-curve
    1233     }
    1234     in.close();
    1235     in.clear();
    1236 
    1237     // read the data for pedestal data
    1238     in_ped.open(input_ped[ch].Data());
    1239     while (in_ped >> da >> dat ){
    1240       hist_scurve_ped[ch].Fill(da,dat); // fill the histogram for pedestal data
    1241     }
    1242     in_ped.close();
    1243     in_ped.clear();
    1244 
    1245     //dac_ped[ch]  = hist_scurve_ped[ch].GetMaximumBin()+dac_min_ped; //absolute value
    1246     dac_ped[ch]  = hist_scurve_ped[ch].GetBinCenter(hist_scurve_ped[ch].GetMaximumBin()); //absolute value
    1247     dac_ped_sum += dac_ped[ch]; //absolute value
    1248     //cerr << "dac_ped[" << ch << "] dac_ped_sum " << dac_ped[ch]
    1249     //<< " " << dac_ped_sum << endl;
    1250    
    1251 
    1252     hist_scurve[ch].Smooth(3); // Smooth the S curve
    1253     Int_t binx=0;
    1254     //hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
    1255     //                                  dac_max-dac_min+1,kHalfMax/2.);
    1256 
    1257     // ??????
    1258     hist_scurve[ch].GetBinWithContent(kHalfMax,binx,kXmin-dac_min+1,
    1259                                       dac_max-dac_min+1,2000);
    1260     dac_half[ch] = hist_scurve[ch].GetBinCenter(binx); //absolute value
    1261     if(dac_half[ch]<kXmin)dac_half[ch]=-1;
    1262    
    1263     // put the histograms of the smoothed S-Curves into two arrays (dac, data)
    1264     //-------------------------------------------------------------------------
    1265     Int_t i=0;
    1266     for(i=0;i<kNdata;i++){
    1267       dac[i]  = hist_scurve[ch].GetBinCenter(i+1);
    1268       data[i] = hist_scurve[ch].GetBinContent(i+1);
    1269     }
    1270    
    1271     /*** compute derive diff on S-Curve ***/
    1272     Double_t diff[kNdata-1], dac_diff[kNdata-1];
    1273     Int_t trigger_diff=0;
    1274     for(i=0;i<kNdata-1;i++){
    1275       diff[i]     = (data[i]-data[i+1])/step_dac;
    1276       dac_diff[i] = dac[i]+step_dac/2;
    1277       hist_diff[ch].Fill(dac_diff[i],diff[i]); // fills differentiated histograms of smoothed s-curves which is the single photoelectron spectrum
    1278       if(dac_diff[i]>=fit_min){
    1279         diff_sum += diff[i]*dac_diff[i];
    1280         //cerr << "diff[" << i << "], dac_diff[" << i << "] "
    1281         //<< diff[i] << " " << dac_diff[i] << endl;
    1282         trigger_diff+=diff[i];
    1283       }
    1284     }
    1285     //mean_diff[ch] = diff_sum/trigger_diff;
    1286     hist_diff[ch].GetXaxis()->SetRangeUser(fit_min,fit_max);
    1287     mean_diff[ch] = hist_diff[ch].GetMean();
    1288     //cerr << "trigger_diff, mean_diff[" << ch << "] "
    1289     //<< trigger_diff << " " << mean_diff[ch] << endl;
    1290    
    1291    
    1292     /******************* GRAPHICS *****************/
    1293 
    1294 
    1295     // Now start the analysis which consist in finding the peak in the differentiated spectrum
    1296     c1->cd(ch+1);   
    1297     /*** fit scurves  ***/
    1298     cerr << "ch " << ch;
    1299 
    1300     //    SDC 24/09/2013 :: Need to declare outside this variable in a compiled program
    1301     TF1 *fit_pol ;
    1302 
    1303 
    1304     // Fit a polynomial curve in a given range on the s-curve directly. This is a kind of DEBUG procedure.
    1305 
    1306     //if(ch==43 || ch==47){ // SDC: 25/09/2013 Why a special channel : the fit pol4 is performed for another DAC range
    1307  /*   if(ch==170){
    1308       fit_pol      = new TF1("fit_pol","pol4",150,fit_max_scurve);
    1309       hist_scurve[ch].Fit(fit_pol,"","",150,fit_max_scurve);
    1310     }else{ */
    1311       fit_pol      = new TF1("fit_pol","pol4",fit_min_scurve,fit_max_scurve);
    1312       hist_scurve[ch].Fit(fit_pol,"","",fit_min_scurve,fit_max_scurve);
    1313    // }
    1314 
    1315    
    1316 
    1317     Double_t fit_par_pol4[5]; // container for fit parameter results
    1318     fit_pol->GetParameters(fit_par_pol4);
    1319     sfit_chiS[ch] =  fit_pol->GetChisquare(); // chi squared
    1320     TF1 *func               = hist_scurve[ch].GetFunction("fit_pol");
    1321     Double_t data_diff_fit;
    1322     for(i=0;i<kNdataFit;i++){
    1323       //if(ch==43 || ch==47){
    1324       if(ch==170){
    1325         dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+150-dac_min;
    1326         data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+150-dac_min); // get the value of the fitted function
    1327       }else{
    1328         dac_fit[i]  = hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min;
    1329         data_fit[i] = func->Eval(hist_scurve[ch].GetBinCenter(i)+fit_min_scurve-dac_min); // get the value of the fitted function
    1330         //cerr << "dac_fit[" << i << "], data_fit[" << i << "] "
    1331         //<< dac_fit[i] << " " << data_fit[i] << endl;
    1332       }
    1333     }
    1334 
    1335     // derive the fitted function on the S-curve to obtain the SPE
    1336     Double_t dac_diff_fit;
    1337     for(i=0;i<kNdataFit-1;i++){
    1338       data_diff_fit = (data_fit[i]-data_fit[i+1])/step_dac;
    1339       hist_diff_fit[ch].Fill(dac_fit[i],data_diff_fit); // Fill differentiated fitted S-curve
    1340     }
    1341 
    1342     TLine *lhalf_max   = new TLine(dac[0],kHalfMax,dac[kNdata-1],kHalfMax);
    1343     TLine *ldac_half   = new TLine(dac_half[ch],0,dac_half[ch],
    1344                                    (Double_t)kDrawMax);
    1345     hist_scurve[ch].SetMaximum(kDrawMax);
    1346     hist_scurve[ch].GetXaxis()->SetRangeUser(kDrawDacMin,kDrawDacMax);
    1347     hist_scurve[ch].Draw();
    1348     lhalf_max->Draw("same"); 
    1349     ldac_half->Draw("same");
    1350     if(dac_half[ch]>0){
    1351       dac_half_sum += dac_half[ch]-dac_ped[ch]; //relative absolute value
    1352       trigger++;
    1353     }
    1354 
    1355     /*** Differential plot from scurve fitting function  ***/
    1356     c5->cd(ch+1);
    1357     hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve,fit_max_scurve-2*step_dac);
    1358     Double_t dac_sfit_max=hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
    1359     //if(dac_sfit_max<=fit_min_scurve+4*step_dac){
    1360     if(dac_sfit_max<=fit_min_scurve+2*step_dac){
    1361       hist_diff_fit[ch].GetXaxis()->SetRangeUser(fit_min_scurve+kFitSadj,fit_max_scurve-2*step_dac);
    1362       dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
    1363       //if(dac_sfit[ch]<=fit_min_scurve+kFitSadj+4*step_dac)dac_sfit[ch]=-1;
    1364       if(dac_sfit[ch]<fit_min_scurve+kFitSadj+step_dac){
    1365         dac_sfit[ch]=-1;
    1366         //dac_sfit[ch]=dac_half[ch];
    1367         flag_fit[ch] = -1;
    1368         cerr << "flag_fit[" << ch << "] " << flag_fit[ch]  << endl;
    1369       }
    1370     }else {
    1371       dac_sfit[ch]  =  hist_diff_fit[ch].GetBinCenter(hist_diff_fit[ch].GetMaximumBin());
    1372     }
    1373 
    1374     if(dac_sfit[ch]>fit_min_scurve+5){
    1375       dac_sfit_sum += dac_sfit[ch]-dac_ped[ch];
    1376       trigger_sfit++;
    1377     }else{
    1378       //dac_sfit[ch] = fit_min_scurve+10;
    1379       dac_sfit[ch] = -1;
    1380     }
    1381     //cerr << "dac_sfit[" << ch << "] " << dac_sfit[ch] << endl;
    1382     TLine *ldac_sfit      = new TLine(dac_sfit[ch],0,dac_sfit[ch],kDrawMaxFit);
    1383     TLine *ldac_sfit_copy = new TLine(dac_sfit[ch],0,dac_sfit[ch],kDrawMax);
    1384     ldac_sfit->SetLineColor(kMagenta);
    1385     ldac_sfit_copy->SetLineColor(kMagenta);
    1386     hist_diff_fit[ch].SetMaximum(kDrawMaxFit);
    1387     hist_diff_fit[ch].SetMinimum(0);
    1388     hist_diff_fit[ch].Draw();
    1389     ldac_sfit->Draw("same");
    1390    
    1391     /**** differential plot ****/
    1392     c2->cd(ch+1);
    1393     //hist_diff[ch].Rebin(8);
    1394     hist_diff[ch].Rebin(3);
    1395     fit_sdiff     = new TF1("fit_sdiff","gaus",dac_min,dac_max);
    1396     //fit_sdiff->SetLineWidth(1);
    1397     //if(ch==8 || ch==18 || ch==44 || ch==56)
    1398     mean_sum += mean_diff[ch]-dac_ped[ch];
    1399     trigger_mean++;
    1400 
    1401     // fit the derivated polynom fit
    1402     Double_t fit_par[3];
    1403     //hist_diff[ch].Fit(fit_sdiff,"RQ0","",fit_min,fit_max);
    1404     hist_diff[ch].Fit(fit_sdiff,"","",fit_min,fit_max);
    1405     hist_diff[ch].GetXaxis()->SetRangeUser(kDrawDacMin,kDrawDacMax);
    1406     hist_diff[ch].SetMinimum(0);
    1407     hist_diff[ch].SetMaximum(kDrawMaxDiff);
    1408     //hist_diff[ch].Rebin();
    1409     hist_diff[ch].Draw();
    1410     fit_sdiff->GetParameters(fit_par); // get the parameters for this fit
    1411     dac_spe[ch] = fit_par[1];
    1412     if(dac_spe[ch]>fit_min){
    1413       dac_spe_sum += dac_spe[ch]-dac_ped[ch];
    1414       //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
    1415       trigger_fit++;
    1416     }else{
    1417       dac_spe[ch] = fit_min+10;
    1418     }
    1419     //cerr << "dac_spe[" << ch << "] " << dac_spe[ch] << endl;
    1420     TLine *ldac_spe    = new TLine(dac_spe[ch],0,dac_spe[ch],
    1421                                  (Double_t)kDrawMaxDiff);
    1422     ldac_spe->SetLineColor(kRed);
    1423     TLine *lmean_diff  = new TLine(mean_diff[ch],0,mean_diff[ch],
    1424                                    (Double_t)kDrawMaxDiff);
    1425     lmean_diff->SetLineColor(kCyan);
    1426     TLine *lfit_min    = new TLine(fit_min,0,fit_min,kDrawMaxDiff);
    1427     lfit_min->SetLineColor(kGray);
    1428 
    1429     ldac_spe->Draw();
    1430     lmean_diff->Draw("same");
    1431     lfit_min->Draw("same");
    1432    
    1433     c3->cd(1);
    1434     if(ch==0){
    1435       hist_scurve[ch].SetMaximum(kDrawMax);
    1436       hist_scurve[ch].Draw();
    1437     }else{
    1438     /*if(ch==8 || ch==18 || ch==44 || ch==56)*/
    1439       hist_scurve[ch].Draw("same");
    1440     }
    1441     if(ch==kDrawCh){
    1442       c3->cd(2);
    1443       hist_scurve[ch].SetMaximum(kDrawMax);
    1444       hist_scurve[kDrawCh].Draw();
    1445       lhalf_max->Draw("same"); 
    1446       ldac_half->Draw("same");
    1447       ldac_sfit_copy->Draw("same");
    1448       c3->cd(3);
    1449       hist_diff[kDrawCh].Draw();     
    1450       ldac_spe->Draw();
    1451       lmean_diff->Draw("same");
    1452       lfit_min->Draw("same");
    1453       c3->cd(4);
    1454       hist_diff_fit[kDrawCh].Draw();
    1455        ldac_sfit->Draw("same");
    1456       }
    1457     /*
    1458     if(ch<32){
    1459       c6->cd(ch/4+1);
    1460       hist_scurve[ch].SetMaximum(kDrawMax);
    1461       if(ch%4==0){
    1462         hist_scurve[ch].Draw();
    1463       }else{
    1464         hist_scurve[ch].Draw("same");
    1465       }
    1466       ldac_sfit_copy->Draw("same");
    1467     }else{
    1468       c7->cd(ch/4-7);
    1469       hist_scurve[ch].SetMaximum(kDrawMax);
    1470       if(ch%4==0){
    1471         hist_scurve[ch].Draw();
    1472       }else{
    1473         hist_scurve[ch].Draw("same");
    1474       }
    1475       ldac_sfit_copy->Draw("same");
    1476       }
    1477     */
    1478    
    1479     c4->cd();
    1480     if(ch==0){
    1481       hist_scurve_ped[ch].SetMaximum(kDrawMaxPed);
    1482       hist_scurve_ped[ch].Draw("l");
    1483     }else{
    1484       hist_scurve_ped[ch].Draw("lsame");
    1485     }
    1486 
    1487 
    1488 if(DEBUG)
    1489 {
    1490     cout << " >>>>>>> END Analysis::ProcessSingleChannel " << endl;     
    1491 }
    1492    
    1493   }// end of process single  channels
     1490        cout << " >>>>>>> END Analysis::ProcessSingleChannel " << endl;
     1491    }
     1492   
     1493}// end of process single  channels
    14941494//---------------------------------------------------------------------------------------------------
    14951495
     
    15001500void Analysis::ComputeDACAverage(Double_t dac_half_spe_fix ,Int_t    kChRef, Bool_t   checkAve )//0:use ref pix, 1:check average
    15011501{
    1502 
    1503 
    1504 if(DEBUG)
    1505 {
    1506   cout << " >>>>>>> START Analysis::ComputeDACAverage " << endl;
    1507  
    1508     cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
    1509     cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
    1510     cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;     
    1511 }
    1512 
    1513  // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
    1514 
    1515  
    1516   mean_diff_ave = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
    1517  
    1518   if(!checkAve){
    1519     dac_ped_ave   = dac_ped[kChRef]; //absolute value
    1520     dac_half_ave  = dac_half[kChRef];
    1521     dac_spe_ave   = dac_spe[kChRef];
    1522     dac_sfit_ave  = dac_sfit[kChRef];
    1523     //Float_t dac_sfit_ave  = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
    1524     mean_diff_ave = mean_diff[kChRef];
    1525   }else{
    1526     dac_ped_ave  = dac_ped_sum/(Float_t)kNch; //absolute value
    1527     dac_half_ave = dac_half_sum/(Float_t)trigger+dac_ped_ave;
    1528     dac_spe_ave  = dac_spe_sum/(Float_t)trigger_fit+dac_ped_ave;
    1529     mean_diff_ave = mean_sum/(Float_t)trigger_mean+dac_ped_ave;
    1530     dac_sfit_ave = dac_sfit_sum/(Float_t)trigger_sfit+dac_ped_ave;
    1531   }
    1532 
    1533 
    1534 // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
    1535 
    1536 
    1537   //trigger_mean=4;
    1538   if(dac_half_spe_fix!=-1) {
    1539     dac_half_spe = dac_half_spe_fix;
    1540   } else {
    1541     //Double_t dac_half_spe = dac_half_ave;
    1542     dac_half_spe = dac_spe_ave;
    1543   }
    1544 
    1545   TLine *ldac_ped_ave   = new TLine(dac_ped_ave,0,dac_ped_ave,
    1546                                    (Double_t)kDrawMax);
    1547   TLine *ldac_ped_ave_copy = new TLine(dac_ped_ave,0,dac_ped_ave,
    1548                                    (Double_t)kDrawMaxDiff);
    1549   TLine *ldac_half_ave  = new TLine(dac_half_ave,0,dac_half_ave,
    1550                                    (Double_t)kDrawMax);
    1551   TLine *ldac_half_spe  = new TLine(dac_half_spe,0,dac_half_spe,
    1552                                    (Double_t)kDrawMax);
    1553   TLine *lmean_diff_ave = new TLine(mean_diff_ave,0,mean_diff_ave,
    1554                                    (Double_t)kDrawMax);
    1555   TLine *ldac_sfit_ave  = new TLine(dac_sfit_ave,0,dac_sfit_ave,
    1556                                    (Double_t)kDrawMax);
    1557   ldac_ped_ave->SetLineColor(kOrange);
    1558   ldac_half_ave->SetLineColor(kRed);
    1559   lmean_diff_ave->SetLineColor(kCyan);
    1560   //ldac_half_spe->SetLineColor(kGreen);
    1561   ldac_sfit_ave->SetLineColor(kMagenta); 
    1562 
    1563   c3->cd(1);
    1564   ldac_ped_ave->Draw("same");
    1565   ldac_half_ave->Draw("same");
    1566   // lhalf_max->Draw("same");  SDC (24/09/2013) :: Not declared
    1567   ldac_half_spe->Draw("same");
    1568   lmean_diff_ave->Draw("same");
    1569   ldac_sfit_ave->Draw("same");
    1570 
    1571   c3->cd(3);
    1572   ldac_ped_ave_copy->SetLineColor(kOrange);
    1573   ldac_ped_ave_copy->Draw("same");
    1574 
    1575 if(DEBUG)
    1576 {
    1577   cout << " >>>>>>> END Analysis::ComputeDACAverage " << endl; 
    1578 }
    1579  
    1580 
     1502   
     1503   
     1504    if(DEBUG)
     1505    {
     1506        cout << " >>>>>>> START Analysis::ComputeDACAverage " << endl;
     1507        
     1508        cout << " \t - dac_half_spe_fix = " << dac_half_spe_fix << endl;
     1509        cout << " \t - kChRef = " << kChRef << " ASIC_C " << endl;
     1510        cout << " \t - checkAve = " << checkAve << " 0:use ref pix, 1:check average " << endl;
     1511    }
     1512   
     1513    // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
     1514   
     1515    
     1516    mean_diff_ave = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
     1517   
     1518    if(!checkAve){
     1519        dac_ped_ave   = dac_ped[kChRef]; //absolute value
     1520        dac_half_ave  = dac_half[kChRef];
     1521        dac_spe_ave   = dac_spe[kChRef];
     1522        dac_sfit_ave  = dac_sfit[kChRef];
     1523        //Float_t dac_sfit_ave  = 252.5; //dac_sfit of ch31 of ASIC_C with the 180 for lower gain pixels
     1524        mean_diff_ave = mean_diff[kChRef];
     1525    }else{
     1526        dac_ped_ave  = dac_ped_sum/(Float_t)kNch; //absolute value
     1527        dac_half_ave = dac_half_sum/(Float_t)trigger+dac_ped_ave;
     1528        dac_spe_ave  = dac_spe_sum/(Float_t)trigger_fit+dac_ped_ave;
     1529        mean_diff_ave = mean_sum/(Float_t)trigger_mean+dac_ped_ave;
     1530        dac_sfit_ave = dac_sfit_sum/(Float_t)trigger_sfit+dac_ped_ave;
     1531    }
     1532   
     1533   
     1534    // SDC 24/09/2013 :: In a compiled program, I have to declare the following variables outside the parenthesis
     1535   
     1536   
     1537    //trigger_mean=4;
     1538    if(dac_half_spe_fix!=-1) {
     1539        dac_half_spe = dac_half_spe_fix;
     1540    } else {
     1541        //Double_t dac_half_spe = dac_half_ave;
     1542        dac_half_spe = dac_spe_ave;
     1543    }
     1544   
     1545    TLine *ldac_ped_ave   = new TLine(dac_ped_ave,0,dac_ped_ave,
     1546                                      (Double_t)kDrawMax);
     1547    TLine *ldac_ped_ave_copy = new TLine(dac_ped_ave,0,dac_ped_ave,
     1548                                         (Double_t)kDrawMaxDiff);
     1549    TLine *ldac_half_ave  = new TLine(dac_half_ave,0,dac_half_ave,
     1550                                      (Double_t)kDrawMax);
     1551    TLine *ldac_half_spe  = new TLine(dac_half_spe,0,dac_half_spe,
     1552                                      (Double_t)kDrawMax);
     1553    TLine *lmean_diff_ave = new TLine(mean_diff_ave,0,mean_diff_ave,
     1554                                      (Double_t)kDrawMax);
     1555    TLine *ldac_sfit_ave  = new TLine(dac_sfit_ave,0,dac_sfit_ave,
     1556                                      (Double_t)kDrawMax);
     1557    ldac_ped_ave->SetLineColor(kOrange);
     1558    ldac_half_ave->SetLineColor(kRed);
     1559    lmean_diff_ave->SetLineColor(kCyan);
     1560    //ldac_half_spe->SetLineColor(kGreen);
     1561    ldac_sfit_ave->SetLineColor(kMagenta);
     1562   
     1563    c3->cd(1);
     1564    ldac_ped_ave->Draw("same");
     1565    ldac_half_ave->Draw("same");
     1566    // lhalf_max->Draw("same");  SDC (24/09/2013) :: Not declared
     1567    ldac_half_spe->Draw("same");
     1568    lmean_diff_ave->Draw("same");
     1569    ldac_sfit_ave->Draw("same");
     1570   
     1571    c3->cd(3);
     1572    ldac_ped_ave_copy->SetLineColor(kOrange);
     1573    ldac_ped_ave_copy->Draw("same");
     1574   
     1575    if(DEBUG)
     1576    {
     1577        cout << " >>>>>>> END Analysis::ComputeDACAverage " << endl;
     1578    }
     1579    
     1580   
    15811581}
    15821582
     
    15861586void Analysis::OutputGain(const TString fname_out)
    15871587{
    1588 
    1589 if(DEBUG)
    1590 {
    1591   cout << " >>>>>>> START Analysis::OutputGain " << endl;
    1592   cout << " \t - fname_out = " << fname_out << endl;   
    1593 }
    1594 
    1595  ofstream out(fname_out.Data());
    1596 
    1597  Int_t    gain_asic_diff, gain_asic_dacHalf, gain_asic_diff_mean, gain_asic_sfit;
    1598   out << "title  ch  g_asic_dacHalf  g_asic_diff  g_asic_diff_mean  g_asic_sfit  dac_half[ch]  dac_ped[ch]  dac_spe[ch]  dac_half_ave dac_half_spe dac_ped_ave  dac_sfit_ave"
    1599        << endl;
    1600 
    1601 
    1602   // loop 3 on all channels to make final computations on the gain
    1603   for(ch=0;ch<kNch;ch++){
    1604     /*
    1605     if(ch==8 || ch==44 || ch==56 || ch==18){
    1606       gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
    1607                                   (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
    1608     }else{
    1609       gain_asic_diff_mean=0;
    1610     }
    1611     */
    1612     if(dac_half[ch]<0){
    1613       gain_asic_dacHalf   = 0;
    1614       gain_asic_diff      = 0;
    1615       gain_asic_diff_mean = 0;
    1616       gain_asic_sfit      = 0;
    1617     }else{
    1618       gain_asic_dacHalf  =(Int_t)((dac_half_ave-dac_ped_ave)/
    1619                                 (dac_half[ch]-dac_ped[ch])*gain_org[ch]);
    1620       gain_asic_diff     =(Int_t)((dac_half_spe-dac_ped_ave)/
    1621                              (dac_spe[ch]-dac_ped[ch])*gain_org[ch]);
    1622       gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
    1623                                   (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
    1624       gain_asic_sfit     =(Int_t)((dac_sfit_ave-dac_ped_ave)/
    1625                              (dac_sfit[ch]-dac_ped[ch])*gain_org[ch]);
    1626       //cerr << ch << " " << gain_org[ch] << endl;
    1627     }
    1628 
    1629    
    1630 
    1631     out << "OUT_TABLE    "
    1632          << ch << "    "
    1633          << gain_asic_dacHalf   << "    "
    1634          << gain_asic_diff      << "    "
    1635          << gain_asic_diff_mean << "    "
    1636          << gain_asic_sfit      << "    "
    1637          << dac_half[ch]        << "    "
    1638          << dac_ped[ch]         << "    "
    1639          << (Int_t)dac_spe[ch]  << "    "
    1640          << dac_half_ave        << "    "
    1641          << dac_half_spe        << "    "
    1642          << dac_ped_ave         << "    "
    1643          << dac_sfit_ave        << "    "
    1644          << sfit_chiS[ch]       << "    "
    1645          << ch                  << "    "
    1646          << flag_fit[ch]       
    1647          << endl;
    1648       std::cout << "channel " << ch << " old gain " << gain_org[ch] << " new gain " << gain_asic_sfit << std::endl;
    1649   }
    1650   out.close();
    1651 
    1652 if(DEBUG)
    1653 {
    1654   cout << " >>>>>>> END Analysis::OutputGain " << endl;
    1655 }
    1656 
    1657 
     1588   
     1589    if(DEBUG)
     1590    {
     1591        cout << " >>>>>>> START Analysis::OutputGain " << endl;
     1592        cout << " \t - fname_out = " << fname_out << endl;
     1593    }
     1594   
     1595    ofstream out(fname_out.Data());
     1596   
     1597    Int_t    gain_asic_diff, gain_asic_dacHalf, gain_asic_diff_mean, gain_asic_sfit;
     1598    out << "title  ch  g_asic_dacHalf  g_asic_diff  g_asic_diff_mean  g_asic_sfit  dac_half[ch]  dac_ped[ch]  dac_spe[ch]  dac_half_ave dac_half_spe dac_ped_ave  dac_sfit_ave"
     1599    << endl;
     1600   
     1601   
     1602    // loop 3 on all channels to make final computations on the gain
     1603    for(ch=0;ch<kNch;ch++){
     1604        /*
     1605         if(ch==8 || ch==44 || ch==56 || ch==18){
     1606         gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
     1607         (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
     1608         }else{
     1609         gain_asic_diff_mean=0;
     1610         }
     1611         */
     1612        if(dac_half[ch]<0){
     1613            gain_asic_dacHalf   = 0;
     1614            gain_asic_diff      = 0;
     1615            gain_asic_diff_mean = 0;
     1616            gain_asic_sfit      = 0;
     1617        }else{
     1618            gain_asic_dacHalf  =(Int_t)((dac_half_ave-dac_ped_ave)/
     1619                                        (dac_half[ch]-dac_ped[ch])*gain_org[ch]);
     1620            gain_asic_diff     =(Int_t)((dac_half_spe-dac_ped_ave)/
     1621                                        (dac_spe[ch]-dac_ped[ch])*gain_org[ch]);
     1622            gain_asic_diff_mean=(Int_t)((mean_diff_ave-dac_ped_ave)/
     1623                                        (mean_diff[ch]-dac_ped[ch])*gain_org[ch]);
     1624            gain_asic_sfit     =(Int_t)((dac_sfit_ave-dac_ped_ave)/
     1625                                        (dac_sfit[ch]-dac_ped[ch])*gain_org[ch]);
     1626            //cerr << ch << " " << gain_org[ch] << endl;
     1627        }
     1628       
     1629       
     1630       
     1631        out << "OUT_TABLE    "
     1632        << ch << "    "
     1633        << gain_asic_dacHalf   << "    "
     1634        << gain_asic_diff      << "    "
     1635        << gain_asic_diff_mean << "    "
     1636        << gain_asic_sfit      << "    "
     1637        << dac_half[ch]        << "    "
     1638        << dac_ped[ch]         << "    "
     1639        << (Int_t)dac_spe[ch]  << "    "
     1640        << dac_half_ave        << "    "
     1641        << dac_half_spe        << "    "
     1642        << dac_ped_ave         << "    "
     1643        << dac_sfit_ave        << "    "
     1644        << sfit_chiS[ch]       << "    "
     1645        << ch                  << "    "
     1646        << flag_fit[ch]       
     1647        << endl;
     1648       
     1649        if(DEBUG>2)
     1650        {
     1651            std::cout << "channel " << ch << " old gain " << gain_org[ch] << " new gain " << gain_asic_sfit << std::endl;
     1652        }
     1653    }
     1654    out.close();
     1655   
     1656    if(DEBUG)
     1657    {
     1658        cout << " >>>>>>> END Analysis::OutputGain " << endl;
     1659       
     1660    }
     1661   
     1662   
    16581663}
    16591664//--------------------------------------------------------------------------------------------------------
     
    16691674void Analysis::Help()
    16701675{
    1671 
    1672   cout << "This is the Help of the PlotScurveAll-diffFit" <<endl;
    1673   cout << " Version 2.0 from September 27 2013 " <<endl;
    1674   cout << " Main author Hiroko Miyamoto " << endl;
    1675   cout << " Adapted for modularity and compilable program by Sylvie Dagoret " << endl;
    1676   cout << " ==> Please from shell execute directly the compiled program :  " <<endl;
    1677   cout << " ==> ../../src/PlotScurvesAllMain  " << endl;
    1678    
     1676   
     1677    cout << "This is the Help of the PlotScurveAll-diffFit" <<endl;
     1678    cout << " Version 2.0 from September 27 2013 " <<endl;
     1679    cout << " Main author Hiroko Miyamoto " << endl;
     1680    cout << " Adapted for modularity and compilable program by Sylvie Dagoret " << endl;
     1681    cout << " ==> Please from shell execute directly the compiled program :  " <<endl;
     1682    cout << " ==> ../../src/PlotScurvesAllMain  " << endl;
     1683    
    16791684}
    16801685//--------------------------------------------------------------------------------------------------------
     
    16841689void Analysis::Menu()
    16851690{
    1686  cout << "This is the Detailed menu information for running the PlotScurveAll-diffFit" <<endl;
    1687 
    1688  cout << " void PlotScurve( " << endl;
    1689  cout << " \t " <<  " Int_t    kDrawCh=31, " << endl;
    1690  cout << " \t " <<  " TString  fname=kFname, " << endl;
    1691  cout << " \t " <<  " TString  fname_ped=kFnamePed, //ASIC_D ped_file " << endl;
    1692  cout << " \t " <<  " Double_t dac_half_spe_fix = -1," << endl;
    1693  cout << " \t " <<  " TString  fname_gain=\"gain-org.txt\"," << endl;
    1694  cout << " \t " <<  " Int_t    fit_min=kFitMin," << endl;
    1695  cout << " \t " <<  " Int_t    fit_max=kFitMax," << endl;
    1696  cout << " \t " <<  " Int_t    fit_min_scurve=kFitMinS," << endl;
    1697  cout << " \t " <<  " Int_t    fit_max_scurve=kFitMaxS," << endl;
    1698  cout << " \t " <<  " Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D " << endl;
    1699  cout << " \t " <<  " Int_t    kChRef=kRefCh,//ASIC_C " << endl;
    1700  cout << " \t " <<  " Bool_t   checkAve=0 //0:use ref pix, 1:check average " << endl;
    1701  cout << " ) " << endl;         
    1702  cout << "Set the parameters on the top of the program depending on your analysis:  " << endl;
    1703  cout << "General policy: parameters starting with a \"k\" at the beginning is fixed gloval parameter and usually set in the header file, never change through the program. " << endl;
    1704  cout << " : " << endl;
    1705  cout << " - kDrawCh: channel just you want to check on display. doens't effect the analysis. " << endl;
    1706  
    1707  cout << " - fname, fname_ped (kFname, kFnamePed): \"YOUR DATA FILE\" without .data"  << endl;
    1708  cout << " See examples in the header file. " << endl;
    1709  cout << " - dac_half_spe_fix: set DAC value only when you want to adjust the gain to particular pulse height (DAC value) instead of either average of the pulse height of reference pixel. Usually just leave it -1. " << endl;
    1710  cout << " - fit_min, fit_max (kFitMin, kFitMax): leave them as default, it is not currently used for the analysis."  << endl;
    1711  cout << " - fit_min_scurve, fit_max_scurve (kFitMinS, kFitMaxS): define the range of s-curve fitting. search proper values by yourself with running root program a few times."  << endl;
    1712  cout << " - kBoard: don't forget to change 0(C) or 1(D) depending on which ASIC you're analysing "  << endl;
    1713  cout << " - kChRef (kRefCh): a pixel whose gain remains 64 thorough the gain adjustment. " << endl;
    1714  cout << " - checkAve: for the first run, set this parameter 1, and find a pixel which remains 64. In the case of parameter 1, all the gains are adjusted to the average pulse height. " << endl;
    1715    cout <<  " Once you choose the reference pixel, set the parameter 0, as well as set kRefCh the pixel number. (better to keep in the header file for the record of each PMTs.) " << endl;
    1716 
    1717    cout << " Actual configuration " << endl;
    1718    cout << " \t " <<  "TString kFname = " << kFname << endl;
    1719    cout << " \t " <<  "TString  kFnamePed = " << kFnamePed  << endl;
    1720    cout << " \t " <<  "Int_t kFitMin = " << kFitMin << endl;
    1721    cout << " \t " <<  "Int_t kFitMax = " << kFitMax << endl;
    1722 
     1691    cout << "This is the Detailed menu information for running the PlotScurveAll-diffFit" <<endl;
     1692   
     1693    cout << " void PlotScurve( " << endl;
     1694    cout << " \t " <<  " Int_t    kDrawCh=31, " << endl;
     1695    cout << " \t " <<  " TString  fname=kFname, " << endl;
     1696    cout << " \t " <<  " TString  fname_ped=kFnamePed, //ASIC_D ped_file " << endl;
     1697    cout << " \t " <<  " Double_t dac_half_spe_fix = -1," << endl;
     1698    cout << " \t " <<  " TString  fname_gain=\"gain-org.txt\"," << endl;
     1699    cout << " \t " <<  " Int_t    fit_min=kFitMin," << endl;
     1700    cout << " \t " <<  " Int_t    fit_max=kFitMax," << endl;
     1701    cout << " \t " <<  " Int_t    fit_min_scurve=kFitMinS," << endl;
     1702    cout << " \t " <<  " Int_t    fit_max_scurve=kFitMaxS," << endl;
     1703    cout << " \t " <<  " Bool_t   kBoard=1, //0:ASIC_C, 1:ASIC_D " << endl;
     1704    cout << " \t " <<  " Int_t    kChRef=kRefCh,//ASIC_C " << endl;
     1705    cout << " \t " <<  " Bool_t   checkAve=0 //0:use ref pix, 1:check average " << endl;
     1706    cout << " ) " << endl;             
     1707    cout << "Set the parameters on the top of the program depending on your analysis:  " << endl;
     1708    cout << "General policy: parameters starting with a \"k\" at the beginning is fixed gloval parameter and usually set in the header file, never change through the program. " << endl;
     1709    cout << " : " << endl;
     1710    cout << " - kDrawCh: channel just you want to check on display. doens't effect the analysis. " << endl;
     1711    
     1712    cout << " - fname, fname_ped (kFname, kFnamePed): \"YOUR DATA FILE\" without .data"  << endl;
     1713    cout << " See examples in the header file. " << endl;
     1714    cout << " - dac_half_spe_fix: set DAC value only when you want to adjust the gain to particular pulse height (DAC value) instead of either average of the pulse height of reference pixel. Usually just leave it -1. " << endl;
     1715    cout << " - fit_min, fit_max (kFitMin, kFitMax): leave them as default, it is not currently used for the analysis."  << endl;
     1716    cout << " - fit_min_scurve, fit_max_scurve (kFitMinS, kFitMaxS): define the range of s-curve fitting. search proper values by yourself with running root program a few times."  << endl;
     1717    cout << " - kBoard: don't forget to change 0(C) or 1(D) depending on which ASIC you're analysing "  << endl;
     1718    cout << " - kChRef (kRefCh): a pixel whose gain remains 64 thorough the gain adjustment. " << endl;
     1719    cout << " - checkAve: for the first run, set this parameter 1, and find a pixel which remains 64. In the case of parameter 1, all the gains are adjusted to the average pulse height. " << endl;
     1720    cout <<  " Once you choose the reference pixel, set the parameter 0, as well as set kRefCh the pixel number. (better to keep in the header file for the record of each PMTs.) " << endl;
     1721   
     1722    cout << " Actual configuration " << endl;
     1723    cout << " \t " <<  "TString kFname = " << kFname << endl;
     1724    cout << " \t " <<  "TString  kFnamePed = " << kFnamePed  << endl;
     1725    cout << " \t " <<  "Int_t kFitMin = " << kFitMin << endl;
     1726    cout << " \t " <<  "Int_t kFitMax = " << kFitMax << endl;
     1727   
    17231728}
    17241729
  • equalization_gain/trunk/src/PlotScurvesAll-diffFitConst.h

    r202 r208  
    5555//TString  kFname="data_s_curves_run94";//T1 ASIC_C it3
    5656//const TString  kFname="data_s_curves_run100";//T1 ASIC_C it3
    57 const TString  kFname="data_s_curves_run104";//T1 ASIC_C it0
     57//const TString  kFname="data_s_curves_run105";//T1 ASIC_C&D it0
     58//const TString  kFname="data_s_curves_run104";//T1 ASIC_C&D it1
     59const TString  kFname="data_s_curves_run105";//T3 ASIC_C&D it1
     60
    5861
    5962/*** pedestal file ***/
     
    7275//TString  kFnamePed="data_s_curves_run94";//T1 ASIC_C it3
    7376//const TString  kFnamePed="data_s_curves_run100";//T1 ASIC_C it3
    74 const TString  kFnamePed="data_s_curves_run104";//T1 ASIC_C it0
     77const TString  kFnamePed="data_s_curves_run105";//T3 ASIC_C&D it1
     78//const TString  kFnamePed="data_s_curves_run104";//T1 ASIC_C&D it1
    7579
    7680// Output file
     
    99103//Int_t kFitSadj = 0;
    100104
    101 const Int_t kRefCh   = 47; //T1 ASIC_D, T1 ASIC_C
     105const Int_t kRefCh   = 4; //T1 ASIC_D, T1 ASIC_C
    102106//Int_t kRefCh   = 32; //(T1 ASIC_C), T2 ASIC_D
    103107//Int_t kRefCh   = 37; //T2 ASIC_C
  • equalization_gain/trunk/src/makeGainTable.sh

    r207 r208  
    1111#export fname=data_s_curves_run91
    1212export fname=data_s_curves_run104
     13#export fname=data_s_curves_run105
    1314#export fname=$@ #standard input
    1415#export fname_o=gain-900V_$fname.txt
     
    2526#export fname_o=gain-1000V-sfit-T1_ASIC_C-th180_th150-ch17-it3.txt
    2627#export fname_o=gain-1000V-sfit-T1_ASIC_C-th200-it3.txt
    27 export fname_o=gain-run_104-T1_ASIC_C-it0.txt
     28export fname_o=gain-run_104-T1_ASIC_D-it1.txt
    2829
    2930
    30 for j in `seq 0 127`; do echo "0"; done > ${fname_o} #ASIC_C only/and ASIC_D
    31 #for j in `seq 0 191`; do echo "0"; done > ${fname_o} #only ASIC_D
     31#for j in `seq 0 127`; do echo "0"; done > ${fname_o} #ASIC_C only/and ASIC_D
     32for j in `seq 0 191`; do echo "0"; done > ${fname_o} #only ASIC_D
    3233
    33 cd ${fname}_C
     34#cd ${fname}_C
     35#grep OUT_TABLE gain-output.dat > gain-output-extract.dat
     36#cat tmp2.dat |awk '{print($3)}' >> ../${fname_o} #dacHalf
     37#cat tmp2.dat |awk '{print($4)}' >> ../${fname_o} #differential
     38#cat gain-output-extract.dat |awk '{print($6)}' >> ../${fname_o} #scurve_fit
     39#cd ../
     40cd ${fname}_D
    3441grep OUT_TABLE gain-output.dat > gain-output-extract.dat
    3542#cat tmp2.dat |awk '{print($3)}' >> ../${fname_o} #dacHalf
    36 #cat tmp2.dat |awk '{print($4)}' >> ../${fname_o} #differential
     43#cat tmp2.dat |awk '{print($5)}' >> ../${fname_o} #differential
    3744cat gain-output-extract.dat |awk '{print($6)}' >> ../${fname_o} #scurve_fit
    3845cd ../
    39 #cd ${fname}_D
    40 #grep OUT_TABLE tmp-v6-th200-it1.data > tmp-v6-th200-it1.dat
    41 #cat tmp2.dat |awk '{print($3)}' >> ../${fname_o} #dacHalf
    42 #cat tmp2.dat |awk '{print($5)}' >> ../${fname_o} #differential
    43 #cat tmp-v6-th200-it1.dat |awk '{print($6)}' >> ../${fname_o} #scurve_fit
    44 #cd ../
    4546
    46 #for j in `seq 0 127`; do echo "0"; done >> ./${fname_o} #ASIC_D only/and ASIC_C
    47 for j in `seq 0 191`; do echo "0"; done >> ./${fname_o} #only ASIC_C
     47for j in `seq 0 127`; do echo "0"; done >> ./${fname_o} #ASIC_D only/and ASIC_C
     48#for j in `seq 0 191`; do echo "0"; done >> ./${fname_o} #only ASIC_C
Note: See TracChangeset for help on using the changeset viewer.