AlignHists.cxx

Go to the documentation of this file.

#include "Alignment/AlignHists.h"

#include "TH2.h"
#include "TFile.h"
#include "TSystem.h"

#include <vector>
#include <iostream>
using std::cerr;
using std::endl;


AlignHists::AlignHists(int nplanes, int nmdlperplane, const char* histname)
    : fIteration(0)
      ,fHistName(histname)
{
    cerr << "AlignHists::AlignHists("
         << nplanes << "," << nmdlperplane << ")\n";

    fNumberOfTracks[0] = new TH2D("NumberOfTracksU",
                                  "Number of tracks per scint module, U view",
                                  nplanes,0,nplanes,
                                  nmdlperplane,0,nmdlperplane);
    fNumberOfTracks[1] = new TH2D("NumberOfTracksV",
                                  "Number of tracks per scint module, V view",
                                  nplanes,0,nplanes,
                                  nmdlperplane,0,nmdlperplane);
    fNumberOfTracksStrips[0] = new TH2D("NumberOfTracksStripsU",
                                  "Number of tracks per strip, U view",
                                  nplanes,0,nplanes,
                                  192,0,192);
    fNumberOfTracksStrips[1] = new TH2D("NumberOfTracksStripsV",
                                  "Number of tracks per strip, V view",
                                  nplanes,0,nplanes,
                                  192,0,192);

    for (int view = 0; view < 2; ++view) {
        char view_char = 'U';
        if (view) view_char = 'V';

        for (int itt = 0; itt < 20; ++ itt) {
            char name[80];
            char title[80];
            sprintf(name,"ResidHist%c%02d",view_char,itt);
            sprintf(title,"Per Module Residuals %c view, iter %02d",
                    view_char,itt);

            fResidHist[view][itt] = new TH2D(name, title,
                                             nplanes,0,nplanes,
                                             nmdlperplane,0,nmdlperplane);
            fResidHist[view][itt]->Sumw2();
            sprintf(name,"ResidSqrHist%c%02d",view_char,itt);
            sprintf(title,"Per Module Residuals**2 %c view, iter %02d",
                    view_char,itt);

            fResidSqrHist[view][itt] = new TH2D(name, title,
                                                nplanes,0,nplanes,
                                                nmdlperplane,0,nmdlperplane);
        }
    }

    cerr << "AlignHists::AlignHists() done\n";
}

AlignHists::AlignHists(const char* filename) 
    : fIteration(0)
      ,fHistName(filename)
{
    cerr << "AlignHists::AlignHists("<<filename<<")\n";

    fNumberOfTracks[0] = new TH2D("NumberOfTracksU",
                                  "Number of tracks per scint module, U view",
                                  500,0,500,8,0,8);
    fNumberOfTracks[1] = new TH2D("NumberOfTracksV",
                                  "Number of tracks per scint module, V view",
                                  500,0,500, 8,0,8);
    fNumberOfTracksStrips[0] = new TH2D("NumberOfTracksStripsU",
                                  "Number of tracks per strip, U view",
                                  500,0,500, 192,0,192);
    fNumberOfTracksStrips[1] = new TH2D("NumberOfTracksStripsV",
                                  "Number of tracks per strip, V view",
                                  500,0,500, 192,0,192);

    for (int view = 0; view < 2; ++view) {
        char view_char = 'U';
        if (view) view_char = 'V';

        for (int itt = 0; itt < 20; ++ itt) {
            char name[80];
            char title[80];
            sprintf(name,"ResidHist%c%02d",view_char,itt);
            sprintf(title,"Per Module Residuals %c view, iter %02d",
                    view_char,itt);

            fResidHist[view][itt] = new TH2D(name, title,500,0,500,8,0,8);
            fResidHist[view][itt]->Sumw2();
            sprintf(name,"ResidSqrHist%c%02d",view_char,itt);
            sprintf(title,"Per Module Residuals**2 %c view, iter %02d",
                    view_char,itt);

            fResidSqrHist[view][itt] = new TH2D(name, title,500,0,500,8,0,8);
        }
    }


}

void AlignHists::LoadFromFile(const char* histname)
{
    TFile* leak = new TFile(histname);
    fNumberOfTracks[0] = (TH2D*)leak->Get("NumberOfTracksU");
    fNumberOfTracks[1] = (TH2D*)leak->Get("NumberOfTracksV");

    for (int view = 0; view < 2; ++view) {
        char view_char = 'U';
        if (view) view_char = 'V';

        for (int itt = 0; itt < 20; ++ itt) {
            char name[80];

            sprintf(name,"ResidHist%c%02d",view_char,itt);
            fResidHist[view][itt] = (TH2D*)leak->Get(name);

            sprintf(name,"ResidSqrHist%c%02d",view_char,itt);
            fResidSqrHist[view][itt] = (TH2D*)leak->Get(name);
        }
    }
    fNumberOfTracksStrips[0] = (TH2D*)leak->Get("NumberOfTracksStripsU");
    fNumberOfTracksStrips[1] = (TH2D*)leak->Get("NumberOfTracksStripsV");    
}

static void add2d(TH2D* h1, TH2D* h2)
{
    int nx = h2->GetNbinsX(), ny = h2->GetNbinsY();
    for (int x=1; x<=nx; ++x) {
        for (int y=1; y<=ny; ++y) {
            h1->SetBinContent(x,y,
                              h1->GetBinContent(x,y) +
                              h2->GetBinContent(x,y));
            h1->SetEntries(h1->GetEntries()+h2->GetEntries());
        }
    }
}

void AlignHists::ReadHists(const char* filename) 
{
    TFile file(filename);

    add2d(fNumberOfTracks[0],(TH2D*)file.Get("NumberOfTracksU"));
    add2d(fNumberOfTracks[1],(TH2D*)file.Get("NumberOfTracksV"));

    add2d(fNumberOfTracksStrips[0],(TH2D*)file.Get("NumberOfTracksStripsU"));
    add2d(fNumberOfTracksStrips[1],(TH2D*)file.Get("NumberOfTracksStripsV"));

    for (int view = 0; view < 2; ++view) {
        char view_char = 'U';
        if (view) view_char = 'V';

        for (int itt = 0; itt < 20; ++ itt) {
            char name[80];

            sprintf(name,"ResidHist%c%02d",view_char,itt);
            add2d(fResidHist[view][itt],(TH2D*)file.Get(name));

            sprintf(name,"ResidSqrHist%c%02d",view_char,itt);
            add2d(fResidSqrHist[view][itt],(TH2D*)file.Get(name));
        }
    }
    
}


AlignHists::~AlignHists()
{
    // dump hists to file?
    if (!fHistName) fHistName = "alignment-hists.root";
    TFile userFile(fHistName,"RECREATE");
    for (int view = 0; view < 2; ++view) {
        for (int itt = 0; itt < 20; ++ itt) {
//            fResidHist[view][itt]->Divide(fNumberOfTracks[view]);
//            fResidSqrHist[view][itt]->Divide(fNumberOfTracks[view]);

            fResidHist[view][itt]->Write();
            fResidSqrHist[view][itt]->Write();
        }
        fNumberOfTracks[view]->Write();
        fNumberOfTracksStrips[view]->Write();
    }
}

void AlignHists::ApplyTrack(PlaneView::PlaneView_t the_view,
                            std::vector<int>& strip,
                            std::vector<int>& plane, std::vector<int>& mdl, 
                            std::vector<double>& resid)
{
    static int count=0;
    ++count;


    int view=0;
    if (the_view == PlaneView::kV) view = 1;
    cerr << "AlignHists::ApplyTrack " << resid.size() << " points,"
         << " iteration=" << fIteration << endl;
    for (unsigned int ind = 0; ind < resid.size(); ++ind) {
        cerr << ".";
        TH2D* hist;

        hist = fResidHist[view][fIteration];
        hist->Fill(plane[ind],mdl[ind],resid[ind]);

        hist = fResidSqrHist[view][fIteration];
        hist->Fill(plane[ind],mdl[ind],resid[ind]*resid[ind]);

        if (!fIteration) {
            fNumberOfTracks[view]->Fill(plane[ind],mdl[ind],1.0);
            fNumberOfTracksStrips[view]->Fill(plane[ind],strip[ind],1.0);
        }
    }
    cerr << endl;
}


ClassImp(AlignHists)

---