Newer
Older
TestStandRepository / Software / ComputeCMN / ComputeCMN.C
//************************************************
// Author: Federica Lionetto
// Created on: 08/01/2014
//************************************************

/*
   ComputeCMN reads a ROOT file (containing the raw data and the data after pedestal subtraction) and creates a ROOT file with the following information:
   - all branches in the Header tree;
   - all branches in the EventInfo tree, one per event (different folders for different CMN subtraction algorithms);
   - float CMN1 (Beetle chip 1), one per event (different folders for different CMN subtraction algorithms);
   - float CMN2 (Beetle chip 2), one per event (different folders for different CMN subtraction algorithms);
   - std::vector<float> ADCCMNSub, one per event (ADC after pedestal and CMN subtraction, different folders for different CMN subtraction algorithms);
   - several distributions, described below.

   The CMN subtraction algorithms are applied after masking the appropriate Beetle channels (cross talk from the header of the Beetle chip, strips and/or Beeetle channels shorted and/or not properly wire bonded).

   Compile with:

   make

   Run with:

   ./ComputeCMN [input ROOT file] [sensor] [additional folder]

where:
- [input ROOT file] is the complete path, including the folder and the filename, of the ROOT file one wants to process;
- [sensor] is the string identifying the sensor (UZHDaughterBoard, Hans410, Hans320, ATLAS, NoSensor);
- [additional folder] is the optional additional folder where the output will be saved.

Example of input ROOT file:
~/LHCb/UTTestBeam/AnalysisResults/CMNStudies/ProcessRawData-filename.root

Example of output ROOT file:
~/LHCb/UTTestBeam/AnalysisResults/CMNStudies/ComputeCMN-filename.root
*/

#include "../Tools/Lib.C"
#include "../Tools/lhcbStyle.C"
#include "../Tools/Style.C"

#include "../Tools/Par.C"

#include "../Tools/Masking.C"

#include "CMNSubMean.C"
#include "CMNSubMean2.C"
#include "CMNSubMeanOddEven.C"
#include "CMNSubMeanxStrips.C"
#include "CMNSubStraightLine.C"

void ComputeCMN(char *filename, char *sensor, char *externalPath=0);

int main(int argc, char *argv[])
{
    getLHCbStyle();
    PersonalStyle();

    if ((argc == 2) && (string(argv[1]) == "--info"))
    {
        cout << "**************************************************" << endl;

        cout << "ComputeCMN reads a ROOT file (containing the raw data and the data after pedestal subtraction) and creates a ROOT file with the following information:" << endl;
        cout << "- all branches in the Header tree;" << endl;
        cout << "- all branches in the EventInfo tree, one per event (different folders for different CMN subtraction algorithms);" << endl;
        cout << "- float CMN1 (Beetle chip 1), one per event (different folders for different CMN subtraction algorithms);" << endl;
        cout << "- float CMN2 (Beetle chip 2), one per event (different folders for different CMN subtraction algorithms);" << endl;
        cout << "- std::vector<float> ADCCMNSub, one per event (ADC after pedestal and CMN subtraction, different folders for different CMN subtraction algorithms);" << endl;
        cout << "- several distributions, described below." << endl;

        cout << " " << endl;

        cout << "The CMN subtraction algorithms are applied after masking the appropriate Beetle channels (cross talk from the header of the Beetle chip, strips and/or Beeetle channels shorted and/or not properly wire bonded)." << endl;

        cout << " " << endl;

        cout << "Compile with:" << endl;

        cout << " " << endl;

        cout << "make" << endl;

        cout << " " << endl;

        cout << "Run with:" << endl;

        cout << " " << endl;

        cout << "./ComputeCMN [input ROOT file] [sensor] [additional folder]" << endl;

        cout << "where:" << endl;
        cout << "- [input ROOT file] is the complete path, including the folder and the filename, of the ROOT file one wants to process;" << endl;
        cout << "- [sensor] is the string identifying the sensor (UZHDaughterBoard, Hans410, Hans320, ATLAS, NoSensor);" << endl;
        cout << "- [additional folder] is the optional additional folder where the output will be saved." << endl;

        cout << "Example of input ROOT file:" << endl;
        cout << "~/LHCb/UTTestBeam/AnalysisResults/CMNStudies/ProcessRawData-filename.root" << endl;

        cout << "Example of output ROOT file:" << endl;
        cout << "~/LHCb/UTTestBeam/AnalysisResults/CMNStudies/ComputeCMN-filename.root" << endl;

        cout << "**************************************************" << endl;

        return 0;
    }
    else if (argc < 3)
    {
        cout << "**************************************************" << endl;

        cout << "Error! Input file missing..." << endl;
        cout << "Please use the following format:" << endl;
        cout << "./ComputeCMN [1] [2] [3]" << endl;
        cout << "with:" << endl;
        cout << "[1] = Input ROOT file, complete path;" << endl;
        cout << "[2] = Sensor (UZHDaughterBoard, Hans410, Hans320, ATLAS, NoSensor);" << endl;
        cout << "[3] = Additional folder, optional." << endl;
        cout << "Type ./ComputeCMN --info for more information." << endl;

        cout << "**************************************************" << endl;

        return 0;  
    }
    else
    {
        cout << "File to process: " << argv[1] << endl;
        cout << "Sensor: " << argv[2] << endl;
        if (argc == 3)
            ComputeCMN(argv[1],argv[2]);
        else if (argc == 4)
            ComputeCMN(argv[1],argv[2],argv[3]);
        else
        {
            cout << "Error! Too many arguments given..." << endl;

            return 0;
        }

        return 0;  
    }
}

void ComputeCMN(char *filename, char *sensor, char *externalPath)
{
    cout << "**************************************************" << endl;
    cout << "Computing CMN..." << endl;
    cout << "**************************************************" << endl;

    // Do not comment this line.
    gROOT->ProcessLine("#include <vector>");

    // Do some fanciness to get the directory right.
    string inputDirectory = (string(filename)).substr(0,(string(filename)).find_last_of('/'));
    // string inputDirectory = "~/LHCb/UTTestBeam/AnalysisResults/CMNStudies";  
    string outputDirectory = inputDirectory;
    // string outputDirectory = "~/LHCb/UTTestBeam/AnalysisResults/CMNStudies";
    string inFilenameNoPathNoExtension = (string(filename)).substr((string(filename)).find_last_of('/')+1);
    inFilenameNoPathNoExtension = inFilenameNoPathNoExtension.substr(0,inFilenameNoPathNoExtension.find_last_of('.'));
    string outFilenameNoPathNoExtension = "ComputeCMN-"+inFilenameNoPathNoExtension.substr(inFilenameNoPathNoExtension.find_first_of('-')+1);
    if (externalPath!=0)
        outputDirectory = string(outputDirectory+"/"+externalPath);
    cout << "The input directory is: " << inputDirectory << endl;
    cout << "The output directory is: " << outputDirectory << endl;
    cout << "The name of the input ROOT file is: " << inFilenameNoPathNoExtension << endl;
    cout << "The name of the output ROOT file is: " << outFilenameNoPathNoExtension << endl;

    // Create the outputDirectory directory if it does not exist.
    string path_to_make = "mkdir -p "+outputDirectory;
    system(path_to_make.c_str());

    cout << "Figures stored in: " << outputDirectory+"/Figures/"+inFilenameNoPathNoExtension << endl;

    // Create a directory named Figures inside the directory named outputDirectory if it does not exist.
    string path_to_make_figures = "mkdir -p "+outputDirectory+"/Figures/"+inFilenameNoPathNoExtension;
    system(path_to_make_figures.c_str());

    TString path_to_figures = (string(path_to_make_figures)).substr((string(path_to_make_figures)).find_last_of(' ')+1);

    string input_ROOT = filename;
    string output_ROOT = outputDirectory+"/"+outFilenameNoPathNoExtension+".root";

    // Open input ROOT file.
    cout << "Open input ROOT file: " << input_ROOT << endl;
    TFile *input = TFile::Open(TString(input_ROOT));

    // Sanity checks for input ROOT file.
    int NEvents;
    TBranch *b_NEvents = 0;

    // Get info from trees.
    TTree *Header = (TTree *)input->Get("Header");
    int EventsHeader = Header->GetEntries();
    if (EventsHeader != 1)
    {
        cout << "Error! The ROOT file has been corrupted..." << endl;

        return ;
    }
    else 
    {
        Header->SetBranchAddress("NEvents",&NEvents,&b_NEvents);

        Header->GetEntry();

        cout << "Number of events: " << NEvents << endl;
    }

    // Get info from trees.
    TTree *EventInfo = (TTree *)input->Get("EventInfo");
    int EventsEventInfo = EventInfo->GetEntries();
    if (EventsEventInfo != NEvents)
    {
        cout << "Error! The ROOT file has been corrupted..." << endl;

        return;
    }

    // Open output ROOT file.
    cout << "Open output ROOT file: " << output_ROOT << endl;
    TFile *output = TFile::Open(TString(output_ROOT),"RECREATE");

    // Copy the Header tree from the input to the output file.
    TTree *HeaderProcessed = Header->CloneTree(0);
    HeaderProcessed->Fill();

    // Copy the histograms with the distribution of the ADC counts of each Beetle channel (TH1F *hADC[N], before pedestal subtraction, and TH1F *hADCPedSub[N], after pedestal subtraction) from the input file to the output file.
    TH1F *hADC[N];
    TH1F *hADCPedSub[N];
    for (int iChannel=0;iChannel<N;++iChannel) {
        hADC[iChannel] = (TH1F *)input->Get(Form("hADC%d",iChannel));
        hADCPedSub[iChannel] = (TH1F *)input->Get(Form("hADCPedSub%d",iChannel));
        
        hADC[iChannel]->Write();
        hADCPedSub[iChannel]->Write();
    }

    // Close input ROOT file.
    input->Close();

    // Masking
    cout << "Masking " << filename << endl;
    Masking(string(sensor));

    // Show the distributions for one particular event.
    int sample = 1000; // Modify it.

    // Configuration of the CMNSubMean CMN subtraction algorithm.
    float th1 = 19.; // Modify it.
    float th2 = 19.; // Modify it.

    CMNSubMean(th1,th2,sample,input_ROOT,output,path_to_figures);

    // Configuration of the CMNSubMean2, CMNSubMeanOddEven, CMNSubMeanxStrips, and CMNSubStraightLine CMN subtraction algorithms.
    float factor = 3.; // Modify it.

    CMNSubMean2(factor,sample,input_ROOT,output,path_to_figures);
    CMNSubMeanOddEven(factor,sample,input_ROOT,output,path_to_figures);

    // Configuration of the CMNSubMeanxStrips CMN subtraction algorithm.
    factor = 3.;
    float x = 64.; // Modify it.

    if (x==0)
        cout << "Error! The number of Beetle channels per group is zero. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else if (NBeetle%(int)x!=0)
        cout << "Error! The number of Beetle channels per group is not a divisor of the total number of Beetle channels. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else 
        CMNSubMeanxStrips(x,factor,sample,input_ROOT,output,path_to_figures);

    factor = 2.2;
    x = 32.; // Modify it.

    if (x==0)
        cout << "Error! The number of Beetle channels per group is zero. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else if (NBeetle%(int)x!=0)
        cout << "Error! The number of Beetle channels per group is not a divisor of the total number of Beetle channels. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else 
        CMNSubMeanxStrips(x,factor,sample,input_ROOT,output,path_to_figures);

    factor = 1.5;
    x = 16.; // Modify it.

    if (x==0)
        cout << "Error! The number of Beetle channels per group is zero. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else if (NBeetle%(int)x!=0)
        cout << "Error! The number of Beetle channels per group is not a divisor of the total number of Beetle channels. The CMNSubMeanxStrips algorithm will not be executed." << endl;
    else 
        CMNSubMeanxStrips(x,factor,sample,input_ROOT,output,path_to_figures);

    factor = 3.;

    CMNSubStraightLine(factor,sample,input_ROOT,output,path_to_figures);

    output->cd();
    output->Write();
    output->Close();

    return;
}