source: trunk/MagicSoft/Mars/macros/status.C@ 2318

/* ======================================================================== *\
!
! *
! * This file is part of MARS, the MAGIC Analysis and Reconstruction
! * Software. It is distributed to you in the hope that it can be a useful
! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
! * It is distributed WITHOUT ANY WARRANTY.
! *
! * Permission to use, copy, modify and distribute this software and its
! * documentation for any purpose is hereby granted without fee,
! * provided that the above copyright notice appear in all copies and
! * that both that copyright notice and this permission notice appear
! * in supporting documentation. It is provided "as is" without express
! * or implied warranty.
! *
!
!
!   Author(s): Thomas Bretz, 4/2003 <mailto:tbretz@astro.uni-wuerzburg.de>
!
!   Copyright: MAGIC Software Development, 2000-2003
!
!
\* ======================================================================== */

/////////////////////////////////////////////////////////////////////////////
//
//  Status - Example how to use the MStatusDisplay
//
//  Use a camera file as input
//
/////////////////////////////////////////////////////////////////////////////

void status(const char *fname="*1947*cosmic*.root")
{
    //
    // Update frequency by default = 1Hz
    //
    MStatusDisplay *d = new MStatusDisplay;

    // Set update time to 5s
    // d->SetUpdateTime(5000);

    // Disable online update
    // d->SetUpdateTime(-1);

    d->SetLogStream(&gLog, kTRUE);            // Disables output to stdout
    gLog.SetOutputFile("status.log", kTRUE);  // Enable output to file
    //gLog.EnableOutputDevice(MLog::eStdout); // Enable output to stdout again

    //
    // Create a empty Parameter List and an empty Task List
    // The tasklist is identified in the eventloop by its name
    //
    MParList  plist;

    MTaskList tlist;
    plist.AddToList(&tlist);

    //
    // The geometry container must be created by yourself to make sure
    // that you don't choose a wrong geometry by mistake
    //
    MGeomCamMagic geomcam;
    plist.AddToList(&geomcam);

    /*
     MSrcPosCam src;
     src.SetXY(1./geomcam.GetConvMm2Deg(), 0);
     plist.AddToList(&src);
     */

    //
    // Now setup the tasks and tasklist:
    // ---------------------------------
    //
    MReadMarsFile read("Events");
    read.DisableAutoScheme();

    // ------------- user change -----------------
    read.AddFile(fname);

    MMcPedestalCopy   pcopy;
    MMcPedestalNSBAdd pnsb;
    MCerPhotCalc      ncalc;

    TArrayS blinds(2);
    blinds[0] = 271;
    blinds[1] = 291;

    MBlindPixelCalc blind;
    //blind.SetPixelIndices(blinds);
    blind.SetUseInterpolation();

    MSigmabarCalc     sgcal;
    MImgCleanStd      clean;
    MHillasCalc       hcalc;
    MHillasSrcCalc    scalc; // !!Preliminary!! Will be removed later!
    MCT1SupercutsCalc calc1;

    MHCamEvent hist("PedestalRms");
    hist.SetType(1);
    plist.AddToList(&hist);

    // -------------------------------------------
    MFillH hfill0("Uncleaned [MHCamEvent]", "MCerPhotEvt");
    MFillH hfill1("Pedestals [MHCamEvent]", "MPedestalCam");
    MFillH hfill2("PedestalRms", "MPedestalCam");
    MFillH hfill3("MHHillas", "MHillas");
    MFillH hfill4("MHHillasExt");
    MFillH hfill5("MHHillasExtSrc [MHHillasExt]", "MHillasSrc");
    MFillH hfill6("MHHillasSrc","MHillasSrc");
    MFillH hfill7("MHNewImagePar","MNewImagePar");
    MFillH hfill8("MHStarMap", "MHillas");
    MFillH hfill9("Cleaned [MHCamEvent]", "MCerPhotEvt");
    MFillH hfill10("MHHadronness", "MHadronness");
    MFillH hfill11("MHSigmaTheta");

    tlist.AddToList(&read);
    tlist.AddToList(&pcopy);
    tlist.AddToList(&pnsb);
    tlist.AddToList(&ncalc);
    tlist.AddToList(&blind);
    tlist.AddToList(&hfill0);
    //tlist.AddToList(&sgcal);
    tlist.AddToList(&clean);
    tlist.AddToList(&hcalc);
    tlist.AddToList(&scalc);
    //tlist.AddToList(&calc1);
    tlist.AddToList(&hfill1);
    tlist.AddToList(&hfill2);
    tlist.AddToList(&hfill3);
    tlist.AddToList(&hfill4);
    tlist.AddToList(&hfill5);
    tlist.AddToList(&hfill6);
    tlist.AddToList(&hfill7);
    tlist.AddToList(&hfill8);
    tlist.AddToList(&hfill9);
    //tlist.AddToList(&hfill10);

    MEvtLoop evtloop;
    evtloop.SetParList(&plist);
    evtloop.SetDisplay(d);

    //
    // Execute your analysis
    //
    if (!evtloop.Eventloop())
        return;

    tlist.PrintStatistics();

    MHCamera &uncl = *((MHCamEvent*)plist.FindObject("Uncleaned"))->GetHistByName();
    MHCamera &hped = *((MHCamEvent*)plist.FindObject("Pedestals"))->GetHistByName();
    MHCamera &hrms = *((MHCamEvent*)plist.FindObject("PedestalRms"))->GetHistByName();

    for (int i=0; i<577; i++)
    {
        if (uncl.IsUsed(i))
        {
            if (uncl[i+1]>uncl.GetMean()+3*uncl.GetRMS())
                cout << "Mean Charge of Pixel-Index #" << i << ": " << uncl[i+1] << " > 3*rms" << endl;
            if (uncl[i+1]==0)
                cout << "Mean Charge of Pixel-Index #" << i << ": " << uncl[i+1] << " <= 0" << endl;
            if (uncl[i+1]<uncl.GetMean()-3*uncl.GetRMS())
                cout << "Mean Charge of Pixel-Index #" << i << ": " << uncl[i+1] << " < 3*rms" << endl;
        }
        if (hped.IsUsed(i))
        {
            if (hped[i+1]>hped.GetMean()+1.5*hped.GetRMS())
                cout << "Mean Pedestal of Pixel-Index #" << i << ": " << hped[i+1] << " > 1.5*rms" << endl;
            if (hped[i+1]==0)
                cout << "Mean Pedestal of Pixel-Index #" << i << ": " << hped[i+1] << " <= 0" << endl;
            if (hped[i+1]<hped.GetMean()-1.5*hped.GetRMS())
                cout << "Mean Pedestal of Pixel-Index #" << i << ": " << hped[i+1] << " < 1.5*rms" << endl;
        }
        if (hrms.IsUsed(i))
        {
            if (hrms[i+1]>hrms.GetMean()+4*hrms.GetRMS())
                cout << "Mean PedestalRMS of Pixel-Index #" << i << ": " << hrms[i+1] << " > 4*rms" << endl;
            if (hrms[i+1]==0)
                cout << "Mean PedestalRMS of Pixel-Index #" << i << ": " << hrms[i+1] << " <= 0" << endl;
            if (hrms[i+1]<hrms.GetMean()-4*hrms.GetRMS())
                cout << "Mean PedestalRMS of Pixel-Index #" << i << ": " << hrms[i+1] << " < 4*rms" << endl;
        }
    }

    //
    // Make sure the display hasn't been deleted by the user while the
    // eventloop was running.
    //
    if ((d = evtloop.GetDisplay()))
    {
        // Save data in a postscriptfile (status.ps)
        //d->SaveAsPS();
        /*
         * ----------- Write status to a root file ------------
         *
         TFile file("status.root", "RECREATE");
         d->Write();
         file.Close();
         delete d;
         */
    }
    /*
     * ----------- Read status from a root file ------------
     *
     TFile file2("status.root", "READ");
     MStatusDisplay *d2 = new MStatusDisplay;
     d2->Read();
     */
}