/*********************************/
/* Compute the hillas parameters */
/*********************************/

Bool_t readDatacards(TString& filename);

// initial and final time slices to be used in signal extraction
const Byte_t hifirst = 0;
const Byte_t hilast  = 13;
const Byte_t lofirst = 3;
const Byte_t lolast  = 12;

// declaration of variables read from datacards
TString  outname;
TString  idirname;
MRunIter pedcaliter;
MRunIter caliter;
MRunIter pediter;
MRunIter datiter;
ULong_t  nmaxevents=999999999;
Short_t  calflag=1;
Float_t  lcore = 3.0;
Float_t  ltail = 1.5;
Int_t islflag = 0;
Float_t  lnew  = 40;
Int_t  kmethod = 1;
Int_t    nfiles = 0;

const TString defaultcard="input.datacard";

//Slow control varialbles
Short_t slowflag=0;
Bool_t  isSlowControlAvailable = kFALSE;
TString hvConfigFile="/mnt/users/jlopez/Mars/Files4Mars/Config/HVSettings_FF35q.conf";
TString continuosLightFile="/local_disk/rootdata/Miscellaneous/Period016/2004_04_16/dc_2004_04_16_04_46_18_22368_Off3c279-2CL100.root";


/*************************************************************/
int makeHillas(TString datacard = defaultcard)
{

  if(!datacard.Length())
    datacard = defaultcard;

  if(!readDatacards(datacard))
    {
      cout << "Error reading datacards. Stoping" << endl;
      return -1;
    }

  // Set the general tasks/containers
  MExtractFixedWindow    extractor;
  extractor.SetRange(hifirst,hilast,lofirst,lolast);

  MCalibrationQECam   qecam;
  MBadPixelsCam       badcam;
  MGeomCamMagic       geomcam;
  MGeomApply          geomapl;

  cout << endl;
  cout << "/****************************************************/" << endl;
  cout << "/* FIRST LOOP: PEDESTAL FOR CALIBRATION COMPUTATION */" << endl;
  cout << "/****************************************************/" << endl;
  cout << endl;

  // If you want to exclude pixels from the beginning, read 
  // an ascii-file with the corr. pixel numbers (see MBadPixelsCam)  
  ifstream fin("badpixels.dat");
  badcam.AsciiRead(fin);
  fin.close();

  MJPedestal pedloop0;
  pedloop0.SetInput(&pedcaliter);
  pedloop0.SetExtractor(&extractor);
  pedloop0.SetBadPixels(badcam);

  if (!pedloop0.Process())
    return;

  cout << endl;
  cout << "/*****************************/" << endl;
  cout << "/* SECOND LOOP: CALIBRATION  */" << endl;
  cout << "/*****************************/" << endl;        
  cout << endl;

  MJCalibration calloop;

  calloop.SetExtractor(&extractor);
  calloop.SetInput(&caliter);
  calloop.SetQECam(qecam);
  calloop.SetBadPixels(pedloop0.GetBadPixels());
  if (!calloop.Process(pedloop0.GetPedestalCam()))
    return;

  cout << endl;
  cout << "/*************************************************/" << endl;
  cout << "/* THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS */" << endl;
  cout << "/*************************************************/" << endl;
  cout << endl;

  // First Compute the pedestals
  MJPedestal pedloop;
  pedloop.SetInput(&pediter);

  if (!pedloop.Process())
    return;

  MPedestalCam pedcam = pedloop.GetPedestalCam();

  // Now convert this pedestals in photons
  MParList  plist3;
  MTaskList tlist3;
  plist3.AddToList(&tlist3);
  
  // containers
  MCerPhotEvt         nphot;
  MPedPhotCam         nphotrms;
  MExtractedSignalCam sigcam;

  plist3.AddToList(&geomcam);
  plist3.AddToList(&pedcam);
  plist3.AddToList(&calloop.GetCalibrationCam());
  plist3.AddToList(&qecam);
  plist3.AddToList(&badcam);
  plist3.AddToList(&sigcam);
  plist3.AddToList(&nphot);
  plist3.AddToList(&nphotrms);

  
  MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault;  
  if(calflag==0)
    calMode=MCalibrate::kNone;
  if(calflag==-1)
    calMode=MCalibrate::kDummy;

  //tasks
  MReadMarsFile read3("Events");
  static_cast<MRead&>(read3).AddFiles(pediter); 
  read3.DisableAutoScheme();
 
//   MReadReports read3;
//   read3.AddTree("Events","MTime.",kFALSE);
//   static_cast<MRead&>(read3).AddFiles(pediter); 

   //  MExtractSignal  extsig;
  //  extsig.SetRange(hifirst,hilast,lofirst,lolast);
  MCalibrate      photcalc(calMode);
  MPedPhotCalc    photrmscalc; 
  
  tlist3.AddToList(&read3);
  tlist3.AddToList(&geomapl);
  tlist3.AddToList(&extractor);
  tlist3.AddToList(&photcalc);
  tlist3.AddToList(&photrmscalc);

  // Create and setup the eventloop
  MEvtLoop evtloop3;
  evtloop3.SetParList(&plist3);
  if (!evtloop3.Eventloop())
    return;
  
  tlist3.PrintStatistics();

  cout << "/**********************************************/" << endl;
  cout << "/* FOURTH LOOP: DATA CALIBRATION INTO PHOTONS */" << endl;
  cout << "/**********************************************/" << endl;

  MParList  plist4;
  MTaskList tlist4;
  plist4.AddToList(&tlist4);
  
  // containers 
  MHillas       hillas;
  MSrcPosCam    source;
  MRawRunHeader runhead;

  MArrivalTimeCam   timecam;
     
  MIslands      isl;
  isl.SetName("MIslands1");
  
  MIslands      isl2;
  isl2.SetName("MIslands2");

  MReportDrive  drive;
  MStarLocalCam starcam;

  if (islflag == 1 || islflag == 2)
    {
      plist4.AddToList(&timecam);
      plist4.AddToList(&isl);
    }
  
  if (islflag == 2)
    {
      plist4.AddToList(&isl2);
    }

  if (isSlowControlAvailable)
    {
      plist4.AddToList(&starcam);
    }
  
  
  plist4.AddToList(&geomcam);
  plist4.AddToList(&pedcam);
  plist4.AddToList(&calloop.GetCalibrationCam());
  plist4.AddToList(&qecam);
  plist4.AddToList(&badcam);
  plist4.AddToList(&nphot);
  plist4.AddToList(&nphotrms);
  plist4.AddToList(&source);
  plist4.AddToList(&hillas);
  plist4.AddToList(&runhead);
  
  //tasks
//   MReadMarsFile read4("Events");
//   static_cast<MRead&>(read4).AddFiles(datiter); 
//   read4.DisableAutoScheme();

  MReadReports read4;
  read4.AddTree("Events","MTime.",kTRUE);
  read4.AddTree("Currents");
  read4.AddTree("Camera");
  read4.AddTree("Drive");
  static_cast<MRead&>(read4).AddFiles(datiter); 

  // Task that use Slow control information 
  MFHVNotNominal fhvnotnominal;
  fhvnotnominal.SetHVNominalValues(hvConfigFile);
  fhvnotnominal.SetMaxNumPixelsDeviated(40);
 
  MContinue hvnotnominal(&fhvnotnominal);

  MCalibrateDC dccalib;
  dccalib.SetFileName(continuosLightFile);

  const Int_t numblind = 5;
  const Short_t w[numblind] = { 47, 124, 470, 475, 571};
  const TArrayS blindpixels(numblind,(Short_t*)w);
  Float_t ringinterest = 100; //[mm]
  Float_t tailcut = 3.5;
  UInt_t integratedevents = 4;

  MFindStars findstars;
  findstars.SetBlindPixels(blindpixels);
  findstars.SetRingInterest(ringinterest);
  findstars.SetDCTailCut(tailcut);
  findstars.SetNumIntegratedEvents(integratedevents);
  findstars.SetMinuitPrintOutLevel(-1);

  MSrcPosFromStars srcposfromstar;
    
  // set bad pixels 
  MBlindPixelCalc   blind;
  MBlindPixelCalc   blind2;
  const Short_t x[16] = {330,395,329,396,389,
                         323,388,322,384,385,
                         386,387,321,320,319,
                         394};
  const TArrayS bp(16,(Short_t*)x);
  blind.SetPixelIndices(bp);
  blind2.SetPixelIndices(bp);
  
  MImgCleanStd      clean(lcore,ltail);

  MArrivalTimeCalc2 timecalc;
  
  MIslandCalc       island;
  island.SetOutputName("MIslands1");

  MIslandClean      islclean(lnew);
  islclean.SetInputName("MIslands1");
  islclean.SetMethod(kmethod);
      
  MIslandCalc       island2;
  island2.SetOutputName("MIslands2");  
  
  
  MHillasCalc       hcalc;
  MHillasSrcCalc    csrc1;
  
  MWriteRootFile write(outname,"RECREATE");
  
  write.AddContainer("MTime"          , "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");

  if (islflag == 1 || islflag == 2)
    write.AddContainer("MIslands1" , "Parameters");
  if (islflag == 2) 
    write.AddContainer("MIslands2" , "Parameters");
  if (isSlowControlAvailable)
    {
      write.AddContainer("MStarLocalCam" , "Parameters");
      write.AddContainer("MReportDrive" , "Parameters");
    }
  
  tlist4.AddToList(&read4);
  tlist4.AddToList(&geomapl);
  if (isSlowControlAvailable)
    {
      tlist4.AddToList(&hvnotnominal,"Events");
      tlist4.AddToList(&dccalib,"Currents");
      tlist4.AddToList(&findstars,"Currents");
      tlist4.AddToList(&srcposfromstar,"Events");
    }
                       
  tlist4.AddToList(&extractor,"Events");
  tlist4.AddToList(&photcalc,"Events");
  //tlist4.AddToList(&blind,"Events");
  tlist4.AddToList(&clean,"Events");

  if (islflag == 1 || islflag == 2)
    {
      tlist4.AddToList(&timecalc,"Events");
      tlist4.AddToList(&island,"Events");
    }

  if (islflag == 2)
    {
      tlist4.AddToList(&islclean,"Events");
      tlist4.AddToList(&island2,"Events");
    }
  
  //tlist4.AddToList(&blind2,"Events");
  tlist4.AddToList(&hcalc,"Events");
  //  tlist4.AddToList(&srcposcalc,"Events");
  tlist4.AddToList(&csrc1,"Events");
  tlist4.AddToList(&write,"Events");

  // Create and setup the eventloop
  MEvtLoop datloop;
  datloop.SetParList(&plist4);

  cout << endl;
  cout << "******************************************************" << endl;
  cout << "* COMPUTING DATA USING EXTRACTED SIGNAL (IN PHOTONS) *" << endl;
  cout << "******************************************************" << endl;
  cout << endl;

  if (!datloop.Eventloop(nmaxevents))
    return;

  tlist4.PrintStatistics();    
}
//-------------------------------------------------------------------------------

Bool_t readDatacards(TString& filename)
{
  ifstream ifun(filename.Data());
  if(!ifun)
    {
      cout << "File " << filename << " not found" << endl;
      return kFALSE;
    }

  TString word;
  
  while(ifun >> word)
    {

      // skip comments
      if(word[0]=='/' && word[1]=='/')
	{
	  while(ifun.get()!='\n'); // skip line
	  continue;
	}

      // number of events
      if(strcmp(word.Data(),"NEVENTS")==0)
	ifun >> nmaxevents;


      // input file directory
      if(strcmp(word.Data(),"IDIR")==0)
	{
	  if(idirname.Length())
	    cout << "readDataCards Warning: overriding input directory file name" << endl;
	  ifun >> idirname;
	}

       // pedestal runs for calibration
      if(strcmp(word.Data(),"PCRUNS")==0)
	{
	  if(pedcaliter.GetNumRuns())
	    cout << "readDataCards Warning: adding pedestal runs for calibration to the existing list" << endl;
	  ifun >> word;
	  pedcaliter.AddRuns(word.Data(),idirname.Data());
	}

      // calibration runs
      if(strcmp(word.Data(),"CRUNS")==0)
	{
	  if(caliter.GetNumRuns())
	    cout << "readDataCards Warning: adding calibration runs to the existing list" << endl;
	  ifun >> word;
	  caliter.AddRuns(word.Data(),idirname.Data());
	}

     // pedestal runs for data
      if(strcmp(word.Data(),"PRUNS")==0)
	{
	  if(pediter.GetNumRuns())
	    cout << "readDataCards Warning: adding pedestal runs for data to the existing list" << endl;
	  ifun >> word;
	  pediter.AddRuns(word.Data(),idirname.Data());
	}

      // data runs
      if(strcmp(word.Data(),"DRUNS")==0)
	{
	  if(datiter.GetNumRuns())
	    cout << "readDataCards Warning: adding data runs to the existing list" << endl;
	  ifun >> word;
	  datiter.AddRuns(word.Data(),idirname.Data());
	}
      
      // output file name
      if(strcmp(word.Data(),"OUTFILE")==0)
	{
	  if(outname.Length())
	    cout << "readDataCards Warning: overriding output file name" << endl;
	  ifun >> outname;
	}

      // calibration flag
      if(strcmp(word.Data(),"CALFLAG")==0)
	ifun >> calflag;

      // cleaning level
      if(strcmp(word.Data(),"CLEANLEVEL")==0)
	{
	  ifun >> lcore;
	  ifun >> ltail;
	}

      if(strcmp(word.Data(),"ISLFLAG")==0)
	{
	  ifun >> islflag;
	}

      // island cleaning 
      if (islflag == 2){
	if(strcmp(word.Data(),"ISLANDCLEAN")==0)
	  {
	    ifun >> kmethod;
	    ifun >> lnew;
	  }
      }

      // slow control data available
      if(strcmp(word.Data(),"SLOWCRT")==0)
	{
	  ifun >> slowflag;
          if (slowflag == 1) isSlowControlAvailable = kTRUE;
	}

      // hv configuration file
      if(strcmp(word.Data(),"HVCONFFILE")==0)
	{
	  ifun >> hvConfigFile;
	}

      // countinous light file to calibrate the dc data
      if(strcmp(word.Data(),"CLFILE")==0)
	{
	  ifun >> continuosLightFile;
	}
    }

   pedcaliter.Reset();
   pediter.Reset();
   caliter.Reset();
   datiter.Reset();
  TString pfile;

  // Dump read values
  cout << "************************************************" << endl;
  cout << "* Datacards read from file " << filename << endl;
  cout << "************************************************" << endl;
  cout << "Pedestal file (s) for calibration: "  << endl;
  while(kTRUE)
  {
      pfile=pedcaliter.Next();
      if (pfile.IsNull())
	  break;
      cout << pfile << endl;
  }
  cout << "Calibration file (s): "  << endl;
  while(kTRUE)
  {
      pfile=caliter.Next();
      if (pfile.IsNull())
	  break;
      cout << pfile << endl;
  }
  cout << "Pedestal file (s) for data: "  << endl;
  while(kTRUE)
  {
      pfile=pediter.Next();
      if (pfile.IsNull())
	  break;
      cout << pfile << endl;
  }
  cout << "Data file (s): "  << endl;
  while(kTRUE)
  {
      pfile=datiter.Next();
      if (pfile.IsNull())
	  break;
      cout << pfile << endl;
  }
  cout << "Maximum number of events: " << nmaxevents << endl;
  cout << "Output file name: " << outname << endl;
  cout << "Calibration flag: " << calflag << endl;
  cout << "Cleaning level: ("<<lcore<<","<<ltail<<")" << endl;
  if (islflag == 1 || islflag == 2)
    cout << "Island calcultation..." << endl;
  if (islflag == 2)
    {
      cout << "Island Cleaning: "<< kmethod <<" method  "<< lnew << " new threshold" << endl;
    }
  cout << "***********" << endl << endl;
  
  if(!pediter.GetNumEntries())
    {
      cout << "No pedestal file name specified" << endl;
      return kFALSE;
    }
  if(!caliter.GetNumEntries())
    {
      cout << "No calibration file name specified" << endl;
      return kFALSE;
    }
  if(!datiter.GetNumEntries())
    {
      cout << "No data file name specified" << endl;
      return kFALSE;
    }
  if(!outname.Length())
    {
      cout << "No output file name specified" << endl;
      return kFALSE;
    }

  if (isSlowControlAvailable)
    {
      cout << "HV Configuratin file: " << hvConfigFile << endl;
      cout << "Continous Light file: " << continuosLightFile << endl;
    }
  

  return kTRUE;
}
