/////////////////////////////////////////////////////////////////////////////
//                                                                         //
// analysis.C                                                              //
//                                                                         //
// Author(s): S.C. Commichau, 6/2004                                       //
//                                                                         //
/////////////////////////////////////////////////////////////////////////////


// This macro calculates source independent and dependent image parameter 
// distributions from data and writes them into a file called hillas.root. 
// To look at the distributions use the macro analysis_read.C
// Up to now no derotation and mispointing correction is taken into account!

const TString inpathcal = "/ihp/MAGIC/Period016/rootdata2/2004_04_22/";
//const TString inpath = "/ihp/MAGIC/Period016/rootdata2/2004_04_23/";
const TString inpath = "/ihp/MAGIC/Period016/rootdata2/2004_04_23/";


const TString dname = "/ihp/MAGIC/Period016/rootdata2/2004_04_23/20040423_2375*_D_OffMrk421-6_E.root";

//const TString dname = "/ihp/MAGIC/Period016/rootdata2/2004_04_23/20040423_2343*_D_Mrk421_E.root";


const Int_t  pedrun  = 23745;//23427;  // Pedestal file
const Int_t  calrun1 = 23205;  // Calibration files
const Int_t  calrun2 = 0;
const Int_t  datrun1 = 23429;  // Data files

const Bool_t usedisplay = kFALSE;

// Range of slices that are taken into account
const Byte_t hifirst = 1;
const Byte_t hilast  = 14;
const Byte_t lofirst = 3;
const Byte_t lolast  = 14;


void analysis()                             
{

    /*******************************************
     * Define some parameters for the analysis *
     *******************************************/
        
    // Change the following line toggle pedestal display 
    Bool_t DrawPedestals = kFALSE;//kTRUE; 
    
    // Set calibration mode
    // kDefault, kNone, kDummy, kFfactor, kPinDiode,...
    MCalibrate::CalibrationMode_t CalMode=MCalibrate::kDefault;

    // Cleaning: Core and Tailcut
    TArrayF clvl(2);
    clvl[0] = 3.5;
    clvl[1] = 3.0;	

    // Set bad pixels
    const Short_t x[16] = {330,395,329,396,389,
			   323,388,322,384,385,
			   386,387,321,320,319,
			   394};

    // Set output filename
    TString oname = "hillas.root";
   
    
    /***************************************
     * Define general tasks and containers *
     ***************************************/

    MCalibrationQECam         qecam;
    MBadPixelsCam             badcam;
    MGeomCamMagic             geomcam;
    MGeomApply                geomapl;
    MExtractFixedWindow       extractor; 
    extractor.SetRange(hifirst, hilast, lofirst, lolast);
    //MExtractFixedWindowPeakSearch extractor;
    //extractor.SetWindows(8,8);
    //MExtractSlidingWindow         extractor;


    // Iterate over runs
    MRunIter pruns;
    MRunIter cruns;
    MRunIter druns;
  
    MReadMarsFile   dread("Events", dname);
    dread.DisableAutoScheme();

    // Pedestal runs
    pruns.AddRun(pedrun,inpath);
    
    // Calibration runs
    if (calrun2==0)
        cruns.AddRun(calrun1,inpathcal);
    else
        cruns.AddRuns(calrun1,calrun2,inpath);

    // Data 
/*    if (1)//datrun2==0)
	druns.AddRun(datrun1,inpath);
    else
	druns.AddRuns(datrun1,datrun2,inpath);//,datrun3,datrun4,datrun5,datrun6,inpath);
*/

    // Look for the pulser colors 
    MCalibrationCam::PulserColor_t color;
 
    if (calrun1 < 20000)
	color = MCalibrationCam::kCT1;
    else
	color = FindColor((MDirIter*)&cruns);
    
    if (color==MCalibrationCam::kNONE)
    {
	TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
	
	while (1)
	{
	    timer.TurnOn();
	    TString input = Getline("Couldn't find the pulser colour, type 'q' to exit, "
				    "'green', 'blue', 'uv' or 'ct1' to go on: ");
	    timer.TurnOff();
	    
	    if (input=="q\n")
		return;
	    
	    if (input=="green\n")
	    {
		color = MCalibrationCam::kGREEN;
		break;
	    }
	    if (input=="blue\n")
	    {
		color = MCalibrationCam::kBLUE;
		break;
	    }
	    if (input=="uv\n" || input=="UV\n")
	    {
		color = MCalibrationCam::kUV;
		break;
	    }
	    if (input=="ct1\n" || input=="CT1\n")
	    {
		color = MCalibrationCam::kCT1;
		break;
	    }
	}
    }
    
 
    /**********************************
     * 1st Loop: Pedestal Computation *
     **********************************/
    
    // If pixels have to be excluded from the beginning, read 
    // an ASCII-file with the corresponding pixel numbers (MBadPixelsCam)
    badcam.AsciiRead("badpixels.txt");
    
    MJPedestal                pedloop;
    pedloop.SetInput(&pruns);   
    pedloop.SetExtractor(&extractor);
    pedloop.SetBadPixels(badcam);
    
    if (!pedloop.Process())
	return;
        
    // Print the content of MPedestalCam
    MPedestalCam *pedcam = (MPedestalCam*)pedloop.GetPedestalCam();
    pedcam->Print(); 
    
    // Draw Pedestals and RMS into a camera - ADC counts
    if (DrawPedestals)
    {
	MHCamera *disp0 = new MHCamera(geomcam, "Pedestals", "");
	MHCamera *disp1 = new MHCamera(geomcam, "PedestalRMS", "");
	
	disp0->SetMinimum(0);
	disp1->SetMinimum(0);
	
	TCanvas *c1 = new TCanvas("Pedestals", "Pedestals", 50, 50, 400, 400);
	TCanvas *c2 = new TCanvas("RMS", "RMS", 460, 50, 400, 400);
	
	c1->SetBorderMode(0);
	c1->SetFillColor(0);
	
	c2->SetBorderMode(0);
	c2->SetFillColor(0);
	
	gStyle->SetOptStat(1101);
	gStyle->SetStatColor(0);
	gStyle->SetStatBorderSize(1);
	
	gPad->SetBorderMode(0);
	
	disp0.SetPrettyPalette();
	disp1.SetPrettyPalette();
	
	c1->cd();
	disp0.Draw();
	
	gPad->cd(1);
	
	disp0.SetCamContent(*pedcam,0);
	disp0.SetCamError(*pedcam,1);
	
	c2->cd();
	disp1.Draw();
	
	disp1.SetCamContent(*pedcam,2);
	disp1.SetCamError(*pedcam,3);
    }
    
   
    /*************************
     * 2nd Loop: Calibration *
     *************************/

    MJCalibration     calloop;
    //calloop.SetColor(color);
    calloop.SetInput(&cruns);
    calloop.SetExtractor(&extractor);
    //calloop.SetQECam(qecam);
    calloop.SetBadPixels(pedloop.GetBadPixels());
    
    // Choose the arrival time Extractor:
    //MExtractTimeHighestIntegral timeext;
    //MExtractTimeSpline timeext;
    // Set Ranges or Windows
    //timeext.SetRange(0,14,6,14);
    //calloop.SetTimeExtractor(&timeext);
    
    // To calibrate times enable the following line
    //Bool_t UseTimes = kFALSE;
    // Bool_t calibrateTimes = kFALSE;
    //calloop.SetRelTimeCalibration(UseTimes);
    
    if (!calloop.Process(pedloop.GetPedestalCam()))
    return;
    
    /**********************************
     * 3rd Loop: Pedestal Calibration *
     **********************************/

    // Create a Parameter and a Task List
    MParList                  plist3;
    MTaskList                 tlist3;
    plist3.AddToList(&tlist3);

    // Parameter Containers
    MCerPhotEvt               evtphot;
    MPedPhotCam               pedphot;
    MExtractedSignalCam       sigcam;
    
    // Add also some containers from previous loops
    plist3.AddToList(&geomcam);
    plist3.AddToList(&pedloop.GetPedestalCam());

    plist3.AddToList(&calloop.GetCalibrationCam());
    plist3.AddToList(&calloop.GetQECam());          // I have to check this! 
    //plist3.AddToList(&calloop.GetRelTimeCam());
    //plist3.AddToList(&calloop.GetBadPixels());
    plist3.AddToList(&badcam);
    plist3.AddToList(&sigcam);
    plist3.AddToList(&evtphot);
    plist3.AddToList(&pedphot);

    
    // Tasks
    MReadMarsFile read3("Events");
    static_cast<MRead&>(read3).AddFiles(pruns);  
    read3.DisableAutoScheme();
    
    //MExtractSignal            extsig;
    //extsig.SetRange(hifirst, hilast, lofirst, lolast);
    MCalibrate                calphot;
    calphot.SetCalibrationMode(CalMode);
    MPedPhotCalc              pedphotcalc;  
    
    tlist3.AddToList(&read3);
    tlist3.AddToList(&geomapl);
    tlist3.AddToList(&extractor);
    //tlist3.AddToList(&extsig);
    tlist3.AddToList(&calphot);
    tlist3.AddToList(&pedphotcalc);
  
    
    // Create and setup the eventloop
    MProgressBar              bar3;
    bar3.SetWindowName("Calibrating Pedestals...");
    
    MEvtLoop                  evtloop3;
    evtloop3.SetProgressBar(&bar3);
    evtloop3.SetParList(&plist3);
    
    if (!evtloop3.Eventloop())
	return;
    
    tlist3.PrintStatistics();


    /****************************
     * 4th Loop: Photonize Data *
     ****************************/

    // Create a Parameter and a Task List
    MParList                  plist4;
    MTaskList                 tlist4;
    plist4.AddToList(&tlist4);

    // Parameter Containers
    MHillas                   hillas;
    //MDCA                      dca;
    MHillasSrc                hillassrc;
    MSrcPosCam                source;
    source.SetX(0.0*189/0.6);
    source.SetY(0.0*189/0.6);
    MRawRunHeader             runhead;
    MArrivalTimeCam           timecam;
    
    // False Source stuff
    MTime                     time;
    MObservatory              obsv;
    MPointingPos              poin;
    
    hillassrc.SetSrcPos(&source);
    

    // Add also some containers from previous loops
    //plist4.AddToList(&timecam);
    plist4.AddToList(&geomcam);
    plist4.AddToList(&pedloop.GetPedestalCam());
    plist4.AddToList(&calloop.GetCalibrationCam());
    //plist4.AddToList(&badcam);
    plist4.AddToList(&calloop.GetQECam());
    //plist4.AddToList(&calloop.GetRelTimeCam());
    //plist4.AddToList(&calloop.GetBadPixels());
    plist4.AddToList(&source);
    plist4.AddToList(&hillas);
    plist4.AddToList(&hillassrc);
    //plist4.AddToList(&dca);    
    plist4.AddToList(&evtphot);
    plist4.AddToList(&pedphot);

    plist4.AddToList(&time);
    plist4.AddToList(&obsv);
    plist4.AddToList(&poin);
    
    
    //MArrivalTime times;
    //plist4.AddToList(&times);
    
    // Tasks
/*  MReadMarsFile read4("Events");
    static_cast<MRead&>(read4).AddFiles(druns);  
    read4.DisableAutoScheme();
*/
    // Set bad pixels
    MBlindPixelCalc           blind;

    const TArrayS bp(16,(Short_t*)x);
    blind.SetPixelIndices(bp);

    MImgCleanStd              clean(clvl[0],clvl[1]);
    MHillasCalc               hilcalc;
    MHillasSrcCalc            hilsrc;
    //MDCACalc                  dcacalc;
    
    
    

    MWriteRootFile write(oname,"RECREATE"); //Alternative options: NEW, UPDATE

    //write.AddContainer("MDCA"         , "Parameters");
    write.AddContainer("MHillas"        , "Parameters");
    write.AddContainer("MHillasSrc"     , "Parameters");
    write.AddContainer("MHillasExt"     , "Parameters");
    write.AddContainer("MNewImagePar"   , "Parameters");
    write.AddContainer("MRawEvtHeader"  , "Parameters");
    write.AddContainer("MRawRunHeader"  , "Parameters");
    write.AddContainer("MConcentration" , "Parameters");
    write.AddContainer("MSrcPosCam"     , "Parameters");

    
    // False Source stuff
    MHFalseSource             fshist;
    fshist.SetAlphaCut(12.5);
    //fshist.SetBgMean(35);
    MFillH                   hfill(&fshist, "MHillas");
    //MFillH                   hfill("MHFalseSource", "MHillas");

    //MFalseSrcCalc     falsecalc;
    //falsecalc.SetAlphaCut(12.5);
    //falsecalc.SetDistCuts(.25,1.0);
    //falsecalc.SetSizeCut(1000);
    //falsecalc.SetWidthCuts(.04,.14);
    //falsecalc.SetLengthCuts(.14,.26);
    //falsecalc.SetOutput("ON.root");//outname);
    //falsecalc.SetRefPoint(0,0);

    //tlist4.AddToList(&read4);
    tlist4.AddToList(&dread);
    tlist4.AddToList(&geomapl);
    tlist4.AddToList(&extractor);
    //tlist4.AddToList(&blind);
    //tlist4.AddToList(&timecalc);
    tlist4.AddToList(&calphot);
    //tlist4.AddToList(&timecal);
    tlist4.AddToList(&clean);
    tlist4.AddToList(&hilcalc);
    //tlist4.AddToList(&dcacalc);
    tlist4.AddToList(&hilsrc);
//    tlist4.AddToList(&falsecalc);
    //tlist4.AddToList(&hfill);
    tlist4.AddToList(&write);  

    // Create and setup the eventloop
    MProgressBar              bar4;
    bar4.SetWindowName("Analyzing Data...");

    MEvtLoop                  evtloop4;
    evtloop4.SetProgressBar(&bar4);
    evtloop4.SetParList(&plist4);

    if (!evtloop4.Eventloop())
	return;

    tlist4.PrintStatistics();

    //fshist.Draw();


}


Bool_t HandleInput()
{
    TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);

    while (1)
    {
        // While reading the input process GUI events asynchronously
        timer.TurnOn();
        TString input = Getline("Type 'q' to exit, <return> to go on: ");
        timer.TurnOff();

        if (input=="q\n")
            return kFALSE;

        if (input=="\n")
            return kTRUE;
    };

    return kFALSE;
}


// Find the correct pulser color, done by M. Gaug
MCalibrationCam::PulserColor_t FindColor(MDirIter* run) 
{
  
  MCalibrationCam::PulserColor_t col = MCalibrationCam::kNONE;

  TString filenames;

  while (!(filenames=run->Next()).IsNull())
    {

      filenames.ToLower();

      if (filenames.Contains("green"))
        if (col == MCalibrationCam::kNONE)
          {
            cout << "Found colour: Green  in " << filenames << endl;
            col = MCalibrationCam::kGREEN;
          }
        else if (col != MCalibrationCam::kGREEN)
          {
            cout << "Different colour found in " << filenames << "... abort" << endl;
            return MCalibrationCam::kNONE;
          }

      if (filenames.Contains("blue"))
        if (col == MCalibrationCam::kNONE)
          {
            cout << "Found colour: Blue  in " << filenames << endl;
            col = MCalibrationCam::kBLUE;
          }
        else if (col != MCalibrationCam::kBLUE)
          {
            cout << "Different colour found in " << filenames << "... abort" << endl;
            return MCalibrationCam::kNONE;
          }

      if (filenames.Contains("uv"))
        if (col == MCalibrationCam::kNONE)
          {
            cout << "Found colour: Uv  in " << filenames << endl;
            col = MCalibrationCam::kUV;
          }
        else if (col != MCalibrationCam::kUV)
          {
            cout << "Different colour found in " << filenames << "... abort" << endl;
            return MCalibrationCam::kNONE;
          }

      if (filenames.Contains("ct1"))
        if (col == MCalibrationCam::kNONE)
          {
            cout << "Found colour: Ct1  in " << filenames << endl;
            col = MCalibrationCam::kCT1;
          }
        else if (col != MCalibrationCam::kCT1)
          {
            cout << "Different colour found in " << filenames << "... abort" << endl;
            return MCalibrationCam::kNONE;
          }
      
    }
  
     
  if (col == MCalibrationCam::kNONE)
    cout <<  "No colour found in filenames of runs: " << ((MRunIter*)run)->GetRunsAsString() 
         << "... abort" << endl;
  
  return col;      
}