//
//-----------------------------------------------------------------------
/// @copyright
///            (c) Copyright 2008 by GATS, Inc.,
///     11864 Canon Blvd, Suite 101, Newport News VA 23606
///
///  All Rights Reserved. No part of this software or publication may be
///  reproduced, stored in a retrieval system, or transmitted, in any form
///  or by any means, electronic, mechanical, photocopying, recording, or
///  otherwise without the prior written permission of GATS, Inc.
///
//-----------------------------------------------------------------------
///
/// @file	Level1Profiles.cpp
///
/// @author	John Burton	<j.c.burton@gats-inc.com>
///
/// @date  	Thu Jan 17 15:25:36 2008
///
//-----------------------------------------------------------------------
//
//-----------------------------------------------------------------------
// Include Files:
//-----------------------------------------------------------------------
//
#include <map>
#include <string>
#include <valarray>
#include <iostream>
#include <sstream>
#include <fstream>

#include "Level1Profiles.h"
#include "GATS_Utilities.hpp"
#include "ChainSawFilter.h"


//
//-----------------------------------------------------------------------
// Defines and Macros:
//-----------------------------------------------------------------------
//


//
//-----------------------------------------------------------------------
// Global Variables:
//-----------------------------------------------------------------------
//

static const double ARCMIN_TO_RADIANS = M_PI / (double)10800.0;
//static EventVar<double> X_Mapping;
//static EventVar<double> Y_Mapping;
//static EventVar<double> X_Mapping(Xname,511.0,-512.0,-1.0);
//static EventVar<double> Y_Mapping(Yname,-511.0,512.0,1.0);
static double det_el = 636.15625;  // now in units of pixels
static double det_az = 518.46875;  // now in units of pixels

static double CSUM_INT    = 347.00; // Center Sum integration time (us)
static double XROI_INT    = 467.67; // X Edge ROI integration time (us)
static double YROI_INT    = 456.83; // Y Edge ROI integration time (us)


//
//-----------------------------------------------------------------------
// Class Methods:
//-----------------------------------------------------------------------
//

void Level1Profiles::initNames(void)
{
  TrackTimeStamp_.setName("SolarTrackerTimes");
  PixelTimeStamp_.setName("SolarImageTimes");
  SunSensorTimeStamp_.setName("SS_ScienceTimeGrid");
  Track_Low_X_.setName("SolarTrackXLow");
  Track_High_X_.setName("SolarTrackXHigh");
  Track_Low_Y_.setName("SolarTrackYLow");
  Track_High_Y_.setName("SolarTrackYHigh");
  
  SolarHighEl_.setName("SolarHighElevationEdge");
  SolarLowEl_.setName("SolarLowElevationEdge");
  SolarHighAz_.setName("SolarHighAzimuthEdge");
  SolarLowAz_.setName("SolarLowAzimuthEdge");
  SolarExtent_.setName("SolarExtent");
  DetectorAz_.setName("DetectorLockdown_AZ");
  DetectorEl_.setName("DetectorLockdown_EL");
  RecordType_.setName("SolarTrackRecordType");
  CenterOffsets_.setName("SolarCenterSumsOffsets");

  NormalizedSolarExtent_.setName("NormalizedSolarExtent");
  
  SumsData_.setName("SolarSums");
  LocsVect_.setName("SolarSumsLocations");
  CenterSumsVect_.setName("SolarCenterSums");
  CenterRowsVect_.setName("SolarCenterRows");
  CenterAvesVect_.setName("CenterAveragePixelValues");
  AvePixVect_.setName("AveragePixelValues");
  CenterOffsetsVect_.setName("SolarCenterSumsOffsetsVect");

}

Level1Profiles::Level1Profiles()
{
}

Level1Profiles::Level1Profiles(SunSensorEvent &sse, Event &L1)
{

  initNames();

  L1.getEventVar("L1_DetectorTimes",   SunSensorTimeStamp_);

  cf_             = sse.getConfigFile();
  section_        = sse.getSection();
  eventNumber_    = sse.getEventNumber();
  srssflag_       = sse.getSRSSFlag();
  RecordType_     = sse.getRecordType(RecordType_);
  TrackTimeStamp_ = sse.getTrackTimeStamp(TrackTimeStamp_);
  PixelTimeStamp_ = sse.getPixelTimeStamp(PixelTimeStamp_);
  Track_Low_X_    = sse.getTrackLoX(Track_Low_X_);
  Track_High_X_   = sse.getTrackHiX(Track_High_X_);
  Track_Low_Y_    = sse.getTrackLoY(Track_Low_Y_);
  Track_High_Y_   = sse.getTrackHiY(Track_High_Y_);
  CenterAvesVect_ = sse.getCenterAvePixVals(CenterAvesVect_);
  CenterSumsVect_ = sse.getCenterSums(CenterSumsVect_);
  CenterRowsVect_ = sse.getCenterSumsRows(CenterRowsVect_);
  CenterOffsetsVect_ = sse.getCenterSumsOffsets(CenterOffsetsVect_);
  LocsVect_          = sse.getAverageLocations(LocsVect_);
  AvePixVect_        = sse.getAveragePixVals(AvePixVect_);

  CSUM_INT = (eventNumber_ > 10732) ? 347.00 : 249.00;
//  Track_Low_X_.Plot("Track_Low_X_","In Level1Profiles");
}


void Level1Profiles::ChainSawFilterAvePixVect(void)
{
  int nelem = AvePixVect_.size();
  int nlen  = AvePixVect_[0].size();
  int minrow = 1023;
  int maxrow = 0;
  int tmp,i,j;
  int istart, istop;
 
  EventVarVect<double> FPASeries;
  EventVarVect<double> AvePixOut;
  EventVar<double> evtmp(nlen);
  EventVar<int>    count(1024);

  std::string cval_name = "Threshold_cval";
  std::string xval_name = "Threshold_xval";
  std::string yval_name = "Threshold_yval";
  double cval_threshold = 0.003;
  double xval_threshold = 0.003;
  double yval_threshold = 0.003;

  if(cf_.ValidEntry(section_,cval_name))
    cval_threshold  = cf_.GetReal(section_,cval_name);
  if(cf_.ValidEntry(section_,xval_name))
    xval_threshold  = cf_.GetReal(section_,xval_name);
  if(cf_.ValidEntry(section_,yval_name))
    yval_threshold  = cf_.GetReal(section_,yval_name);

  evtmp.setValueZero();
  for(i=0;i<nelem;i++)
    AvePixOut.addEventVar(evtmp);

  for(i=0;i<1024; i++)
    FPASeries.addEventVar(evtmp);

//
// Do X Edges
//
  istart = 0;
  istop  = 14;
  count.setValueZero();
  for(i=0;i<1024;i++)
    FPASeries[i].setValueToMissing();
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      tmp = LocsVect_[i][j];
      FPASeries[tmp][j] = AvePixVect_[i][j];
      count[tmp] = count[tmp] + 1;
    }
  }
  int maxcount = count.max();
  for(i=0;i<1024;i++)
  {
    if(count[i] > 2)
    {
      FPASeries[i] = applyChainSawFilter(FPASeries[i],xval_threshold);
    }
  }
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      tmp = LocsVect_[i][j];
      AvePixVect_[i][j] = FPASeries[tmp][j];
    }
  }
//
// Do y Edges
//
  istart = 14;
  istop  = 24;
  count.setValueZero();
  for(i=0;i<1024;i++)
    FPASeries[i].setValueToMissing();
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      tmp = LocsVect_[i][j];
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      FPASeries[tmp][j] = AvePixVect_[i][j];
      count[tmp] = count[tmp] + 1;
    }
  }
  for(i=0;i<1024;i++)
  {
    if(count[i] > 2)
    {
      FPASeries[i] = applyChainSawFilter(FPASeries[i],xval_threshold);
    }
  }
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      tmp = LocsVect_[i][j];
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      AvePixVect_[i][j] = FPASeries[tmp][j];
    }
  }
//
// Do Center Sums
//
  istart = 24;
  istop  = 31;
  count.setValueZero();
  for(i=0;i<1024;i++)
    FPASeries[i].setValueToMissing();
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      tmp = LocsVect_[i][j];
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      FPASeries[tmp][j] = AvePixVect_[i][j];
      count[tmp] = count[tmp] + 1;
    }
  }
  for(i=0;i<1024;i++)
  {
    if(count[i] > 2)
    {
      FPASeries[i] = applyChainSawFilter(FPASeries[i],xval_threshold);
    }
  }
  for(i=istart;i<istop;i++)
  {
    for(j=0;j<nlen;j++)
    {
      tmp = LocsVect_[i][j];
      if(tmp < 0) 
	tmp = 0;
      else if (tmp > 1023)
	tmp = 1023;
      AvePixVect_[i][j] = FPASeries[tmp][j];
    }
  }
}



void Level1Profiles::doChainSawFilter(void)
{
  EventVarVect<double> evv;
  std::string lox_name = "Threshold_lox";
  std::string hix_name = "Threshold_hix";
  std::string loy_name = "Threshold_loy";
  std::string hiy_name = "Threshold_hiy";
  double lox_threshold = 0.000015; 
  double hix_threshold = 0.000015; 
  double loy_threshold = 0.0002; 
  double hiy_threshold = 0.0002; 

  if(cf_.ValidEntry(section_,lox_name))
    lox_threshold  = cf_.GetReal(section_,lox_name);
  if(cf_.ValidEntry(section_,hix_name))
    hix_threshold  = cf_.GetReal(section_,hix_name);
  if(cf_.ValidEntry(section_,loy_name))
    loy_threshold  = cf_.GetReal(section_,loy_name);
  if(cf_.ValidEntry(section_,hiy_name))
    hiy_threshold  = cf_.GetReal(section_,hiy_name);

  evv.addEventVar(Track_Low_X_);
  std::cerr << "doChainSawFilter::Track_Low_X\n";
  Track_Low_X_ = applyChainSawFilter(Track_Low_X_,lox_threshold);
  evv.addEventVar(Track_Low_X_);
  if(cf_.GetFlag(section_,"PlotTrackingData"))
    evv.Plot("Low_X_Track_","Level1Profiles::doChainSawFilter");
  evv.clear();
  evv.addEventVar(Track_High_X_);
//  std::cerr << "doChainSawFilter::Track_High_X\n";
  Track_High_X_ = applyChainSawFilter(Track_High_X_,hix_threshold);
  evv.addEventVar(Track_High_X_);
  if(cf_.GetFlag(section_,"PlotTrackingData"))
    evv.Plot("High_X_Track_","Level1Profiles::doChainSawFilter");
  evv.clear();
//  std::cerr << "doChainSawFilter::Track_Low_Y\n";
  Track_Low_Y_ = applyChainSawFilter(Track_Low_Y_,loy_threshold);
//  std::cerr << "doChainSawFilter::Track_High_Y\n";
  Track_High_Y_ = applyChainSawFilter(Track_High_Y_,hiy_threshold);
//  Track_Low_X_.Plot("Track_Low_X_","After Chainsaw Filter");
  if(cf_.GetFlag(section_,"PlotTrackingData"))
  {
    EventVarVect<double> evv;
    evv.addEventVar(Track_Low_X_);
    evv.addEventVar(Track_High_X_);
    evv.addEventVar(Track_Low_Y_);
    evv.addEventVar(Track_High_Y_);
    evv.Plot("Filtered Tracking Data","Level1Profiles::doChainSawFilter");
    evv.clear();
  }


}



void Level1Profiles::transformCoordinates(Event &SD)
{

  std::string fname;

  std::string Xname  = "X_Pixel_Spacing";
  std::string Xfile  = "X_MappingFile";

  std::string Yname  = "Y_Pixel_Spacing";
  std::string Yfile  = "Y_MappingFile";

  std::string det_az_name;
  std::string det_el_name;

  if(srssflag_ == 'r')
  {
    det_az_name = "Detector_lockdown_az_SR";
    det_el_name = "Detector_lockdown_el_SR";
  }
  else if(srssflag_ == 's')
  {
    det_az_name = "Detector_lockdown_az_SS";
    det_el_name = "Detector_lockdown_el_SS";
  }
  else
  {
    det_az_name = "Detector_lockdown_az";
    det_el_name = "Detector_lockdown_el";
  }
  EventVar<double> X_Mapping(Xname,511.0,-512.0,-1.0);
  EventVar<double> Y_Mapping(Yname,-511.0,512.0,1.0);

  if(SD.EventVarExists(Xname))
    SD.getEventVar(Xname,X_Mapping);
  else
  {
    if(cf_.ValidEntry(section_,Xfile))
    {
      fname  = cf_.GetStr(section_,Xfile);
      std::cerr << "Reading " << fname << std::endl;
      X_Mapping.readVar(fname);
    }
    X_Mapping = X_Mapping * ARCMIN_TO_RADIANS;
    SD.addEventVar(X_Mapping);
  }

  if(SD.EventVarExists(Yname))
    SD.getEventVar(Yname,Y_Mapping);
  else
  {
    if(cf_.ValidEntry(section_,Yfile))
    {
      fname  = cf_.GetStr(section_,Yfile);
      std::cerr << "Reading " << fname << std::endl;
      Y_Mapping.readVar(fname);
    }
    Y_Mapping = Y_Mapping * ARCMIN_TO_RADIANS;
    SD.addEventVar(Y_Mapping);
  }

  if(cf_.ValidEntry(section_,det_az_name))
    det_az  = cf_.GetReal(section_,det_az_name);

  if(cf_.ValidEntry(section_,det_el_name))
    det_el  = cf_.GetReal(section_,det_el_name);

  det_az = Y_Mapping.Interpol_ndx(det_az);
  det_el = X_Mapping.Interpol_ndx(det_el);
  SolarHighEl_ = X_Mapping.VInterpol_ndx(Track_Low_X_);
  SolarLowEl_  = X_Mapping.VInterpol_ndx(Track_High_X_);
  SolarHighAz_ = Y_Mapping.VInterpol_ndx(Track_High_Y_);
  SolarLowAz_  = Y_Mapping.VInterpol_ndx(Track_Low_Y_);
  for (int i=0; i<CenterOffsetsVect_.size(); i++)
    CenterOffsetsVect_[i] = X_Mapping.VInterpol_ndx(CenterOffsetsVect_[i]);
//  SolarHighEl_.Plot("SolarHighEl_","After TransformCoordinates");
}



void Level1Profiles::determineDetectorLockDownPosition(void)
{
  EventVar<double> SunCenterAz;

  SunCenterAz = (SolarLowAz_ + SolarHighAz_) / 2.0;
  DetectorAz_ = SunCenterAz  - det_az;
  DetectorEl_ = SolarHighEl_ - det_el;
}


void Level1Profiles::determineCenterSumsElevation(void)
{

  for(int i = 0; i<CenterOffsetsVect_.size(); i++)
      CenterOffsetsVect_[i] = det_el - CenterOffsetsVect_[i];
}


void Level1Profiles::determineSolarExtent(void)
{
  double exoExtent;
//  double stddev;
  std::string str;
  SolarExtent_ = SolarHighEl_ - SolarLowEl_;

  valarray<bool> good  = ((TrackTimeStamp_ >= 60.0) &&
			  (TrackTimeStamp_ <= 80.0));
  
//  SolarExtent_.Plot("Solar Extent","FPACorrection");
  EventVar<double> tmp = SolarExtent_;
//  EventVar<double> mstd;
//  mstd = tmp[std::slice(1700,500,1)];
//  double mean   = mstd.Mean();
//  double stddev = mstd.Stddev();
//  std::cerr << "SolarExtent  mean, std dev: ";
//  std::cerr << mean << ", " << stddev << "\n"; 

  tmp = tmp[good];

//  SolarExtent_.Plot("Solar Extent","");
  exoExtent = tmp.max();

//  str = "Event " + GATS_Utilities::ConvertToString(eventNumber_) + 
//        " SD = " + GATS_Utilities::ConvertToString(tmp.Stddev());
//  std::cerr << "tmp.stddev() = " << tmp.Stddev() << "\n";
  NormalizedSolarExtent_ = SolarExtent_ / exoExtent;
}


EventVar<double> Level1Profiles::getSolarExtent(void)
{
  EventVar<double> SolarAzExtent;

  SolarAzExtent = (SolarHighAz_ + SolarLowAz_) / 2.0;
//  return SolarExtent_;
  return SolarAzExtent;
}


void Level1Profiles::adjustToScienceTimeStamp(Event& L1)
{

  std::string offset_name = "SunSensorTimeShift";
  std::string str, tmpstr;
  std::stringstream oss;
  double SS_TimeShift = 0.0;

  if(cf_.ValidEntry(section_,offset_name))
    SS_TimeShift  = cf_.GetReal(section_,offset_name);
 
  oss << eventNumber_;
  str = "Event " + oss.str();

  EventVar<double> xloc_ts, csum_ts, hiy_ts, loy_ts, xsum_ts;

  xloc_ts = TrackTimeStamp_ - SS_TimeShift;
  csum_ts = TrackTimeStamp_ - SS_TimeShift;
  hiy_ts  = TrackTimeStamp_ - SS_TimeShift;
  loy_ts  = TrackTimeStamp_ - SS_TimeShift;
  xsum_ts = TrackTimeStamp_ - SS_TimeShift;

//  xloc_ts = TrackTimeStamp_ - 0.05;
//  csum_ts = TrackTimeStamp_ - 0.04;
//  hiy_ts  = TrackTimeStamp_ - 0.08;
//  loy_ts  = TrackTimeStamp_ - 0.07;
//  xsum_ts = TrackTimeStamp_ - 0.02;

  SolarHighEl_ = SolarHighEl_.VInterpol(xloc_ts, SunSensorTimeStamp_);
  SolarLowEl_  = SolarLowEl_.VInterpol (xloc_ts, SunSensorTimeStamp_);
  SolarHighAz_ = SolarHighAz_.VInterpol(hiy_ts, SunSensorTimeStamp_);
  SolarLowAz_  = SolarLowAz_.VInterpol (loy_ts, SunSensorTimeStamp_);
  SolarExtent_ = SolarExtent_.VInterpol(xloc_ts, SunSensorTimeStamp_);
  NormalizedSolarExtent_ = NormalizedSolarExtent_.VInterpol(xloc_ts, SunSensorTimeStamp_);

  DetectorAz_  = DetectorAz_.VInterpol (xloc_ts, SunSensorTimeStamp_);
  DetectorEl_  = DetectorEl_.VInterpol (xloc_ts, SunSensorTimeStamp_);

  CenterAvesVect_ = CenterAvesVect_.VInterpol(csum_ts,SunSensorTimeStamp_);

//  SolarHighEl_.Plot("SolarHighEl_","After AdjustToScienceTimeGrid");
}


void Level1Profiles::copyDataToEvent(Event& eventobj)
{
  initNames();
  eventobj.addEventVar(TrackTimeStamp_);
  eventobj.addEventVar(Track_Low_X_);
  eventobj.addEventVar(Track_High_X_);
  eventobj.addEventVar(Track_Low_Y_);
  eventobj.addEventVar(Track_High_Y_);
  eventobj.addEventVar(SolarHighEl_);
  eventobj.addEventVar(SolarLowEl_);
  eventobj.addEventVar(SolarHighAz_);
  eventobj.addEventVar(SolarLowAz_);
  eventobj.addEventVar(SolarExtent_);
  eventobj.addEventVar(DetectorEl_);
  eventobj.addEventVar(DetectorAz_);
  eventobj.addEventVar(PixelTimeStamp_);
  eventobj.addEventVar(SumsData_);
  eventobj.addEventVar(NormalizedSolarExtent_);
  eventobj.addEventVar(CenterOffsets_);
  eventobj.addEventVar(CenterAvesVect_);
//  eventObj.addEventVar(CenterOffsetsVect_);
  eventobj.addEventVar(CenterSumsVect_);
  eventobj.addEventVar(SunSensorTimeStamp_);
  eventobj.addEventVar(LocsVect_);
  eventobj.addEventVar(AvePixVect_);
}


void Level1Profiles::correctCenterSumsLocations()
{
  std::stringstream oss;
  std::string str, tmpstr;
  std::map<int, EventVar<double> > CenterAveMap;
  std::map<int, EventVar<double> > CenterOffsetMap;
  std::map<int, EventVar<double> >::iterator iter;
  std::vector<int> rowvect;
  int nvars = CenterAvesVect_.size();
  int len   = TrackTimeStamp_.size();

  valarray<bool> good  = ((TrackTimeStamp_ >= 70.0) &&
			  (TrackTimeStamp_ <= 71.0));

  EventVar<int> Rows;
  EventVar<double> Aves;
  EventVar<double> Offs;
  EventVar<double> tmp(len);

  int row;
  CenterOffsets_.resize(nvars);
  for(int j=0; j<nvars; j++)
  {
    Rows = CenterRowsVect_[j];
    Rows = Rows[good];
    row  = Rows.max();
    rowvect.push_back(row);
    str = "CenterAve_" + GATS_Utilities::ConvertToString(row);
    tmp.setName(str);
    tmp.setValueToMissing();
    CenterAveMap[row] = tmp;
    str = "CenterOffset_" + GATS_Utilities::ConvertToString(row);
    tmp.setName(str);
    tmp.setValueToMissing();
    CenterOffsetMap[row] = tmp;
  }
  for(int j=0; j<nvars; j++)
  {
    Rows = CenterRowsVect_[j];
    Aves = CenterAvesVect_[j];
    Offs = CenterOffsetsVect_[j][good];
    CenterOffsets_[j] = Offs.Mean();
    for(int i=0; i<Rows.size(); i++)
    {
      row = Rows[i];
      iter = CenterAveMap.find(row);
      if(iter != CenterAveMap.end())
      {
	CenterAveMap[row][i]    = Aves[i];
      }
    }
  }
  CenterAvesVect_.clear();
  CenterAvesVect_.setName("AveragePixelValues");
  CenterOffsets_.setName("SolarCenterSumsOffsets");
  for(int j=0; j<rowvect.size(); j++)
  {
    row  = rowvect[j];
    Aves = CenterAveMap[row];
//    Offs = (*CenterOffsetMap[row]);
    CenterAvesVect_.addEventVar(Aves);
//    CenterOffsetsVect_->addEventVar(Offs);
  }
  
}


void Level1Profiles::correctSumsLocations()
{
  std::stringstream oss;
  std::string str, tmpstr;
  std::map<int, EventVar<double> > AveMap;
  std::map<int, EventVar<double> > LocationMap;
  std::map<int, EventVar<double> >::iterator iter;
  std::vector<int> rowvect;
  int nvars = AvePixVect_.size();
//  int len   = sstVec_.size();
  int len   = TrackTimeStamp_.size();

  valarray<bool> good  = ((TrackTimeStamp_ >= 70.0) &&
			  (TrackTimeStamp_ <= 71.0));

  EventVar<int> Rows;
  EventVar<double> Aves;
  EventVar<double> Offs;
  EventVar<double> tmp(len);

  int row;
  for(int j=0; j<nvars; j++)
  {
    Rows = LocsVect_[j];
    int rowmax = Rows.max();
    int rowmin = Rows.min();
    std::cerr << "Row[" << j << "] = " << rowmax << ", " << rowmin << " :: ";
    Rows = Rows[good];
    rowmax = Rows.max();
    rowmin = Rows.min();
    std::cerr << rowmax << ", " << rowmin << std::endl;
    row  = Rows.max();
    rowvect.push_back(row);
    str = "Ave_" + GATS_Utilities::ConvertToString(row);
    tmp.setName(str);
    tmp.setValueToMissing();
    AveMap[row] = tmp;
    str = "Location_" + GATS_Utilities::ConvertToString(row);
    tmp.setName(str);
    tmp.setValueToMissing();
    LocationMap[row] = tmp;
  }
  for(int j=0; j<nvars; j++)
  {
    Rows = LocsVect_[j];
    Aves = AvePixVect_[j];
    for(int i=0; i<Rows.size(); i++)
    {
      row = Rows[i];
      iter = AveMap.find(row);
      if(iter != AveMap.end())
      {
	AveMap[row][i]    = Aves[i];
      }
    }
  }
  AvePixVect_.clear();
  AvePixVect_.setName("AveragePixelValues");
  for(int j=0; j<rowvect.size(); j++)
  {
    row  = rowvect[j];
    Aves = AveMap[row];
    AvePixVect_.addEventVar(Aves);
  }
 
  oss << eventNumber_;
  str = "Event " + oss.str();
//  AvePixVect_->Plot("Average PixelValues",str,3,4);
//  CenterOffsetsVect_->Plot("Center Sums Offsets",str);
}


void Level1Profiles::smoothData(void)
{
  int len = TrackTimeStamp_.size();
  int i,j,n;
  int maf_points = 1;
  int maf_passes = 0;
  std::vector<std::string> maf_variables;
  std::string MAFVars   = "MAF_Variables";
  std::string MAFPoints = "MAF_Points";
  std::string MAFPasses = "MAF_Passes";
  if(cf_.ValidEntry(section_,MAFPoints))
  {
    maf_points  = cf_.GetInt(section_,MAFPoints);
  }

  if(cf_.ValidEntry(section_,MAFPasses))
  {
    maf_passes  = cf_.GetInt(section_,MAFPasses);
  }

//  if(cf_.ValidEntry(section_,MAFVars))
//  {
//    maf_variables = cf_.GetStrValues(section_,MAFVars);
//    for(n = 0; n < maf_variables.size(); n++)
//    {
//      std::string str = maf_variables[n];
//      switch(VarMap[str])
//      {
//      case evLowX:
//	for(i=0; i<maf_passes; i++)
//	  (*Track_Low_X_) = Track_Low_X_->Smooth(maf_points);
//	Track_Low_X_->setName(str);
//	break;
//      case evHighX:
//	for(i=0; i<maf_passes; i++)
//	  (*Track_High_X_) = Track_High_X_->Smooth(maf_points);
//	Track_High_X_->setName(str);
//	break;
//      case evLowY:
//	for(i=0; i<maf_passes; i++)
//	  (*Track_Low_Y_) = Track_Low_Y_->Smooth(maf_points);
//	Track_Low_Y_->setName(str);
//	break;
//      case evHighY:
//	for(i=0; i<maf_passes; i++)
//	  (*Track_High_Y_) = Track_High_Y_->Smooth(maf_points);
//	Track_High_Y_->setName(str);
//	break;
//      case evHighEl:
//	for(i=0; i<maf_passes; i++)
//	  (*SolarHighEl_) = SolarHighEl_->Smooth(maf_points);
//	SolarHighEl_->setName(str);
//	break;
//      case evLowEl:
//	for(i=0; i<maf_passes; i++)
//	  (*SolarLowEl_) = SolarLowEl_->Smooth(maf_points);
//	SolarLowEl_->setName(str);
//	break;
//      case evHighAz:
//	for(i=0; i<maf_passes; i++)
//	  (*SolarHighAz_) = SolarHighAz_->Smooth(maf_points);
//	SolarHighAz_->setName(str);
//	break;
//      case evLowAz:
//	for(i=0; i<maf_passes; i++)
//	  (*SolarLowAz_) = SolarLowAz_->Smooth(maf_points);
//	SolarLowAz_->setName(str);
//	break;
//      case evExtent:
//	for(i=0; i<maf_passes; i++)
//	  (*SolarExtent_) = SolarExtent_->Smooth(maf_points);
//	SolarExtent_->setName(str);
//	break;
//      case evDetAz:
//	for(i=0; i<maf_passes; i++)
//	  (*DetectorAz_) = DetectorAz_->Smooth(maf_points);
//	DetectorAz_->setName(str);
//	break;
//      case evDetEl:
//	for(i=0; i<maf_passes; i++)
//	  (*DetectorEl_) = DetectorEl_->Smooth(maf_points);
//	DetectorEl_->setName(str);
//	break;
//      case evSums:
//	for(j=0; j<SumsData_->size(); j++)
//	{
//	  for(i=0; i<maf_passes; i++)
//	    (*SumsData_)[j] = (*SumsData_)[j].Smooth(maf_points);
//	  std::string tmpstr = str + GATS_Utilities::ConvertToString(j);
//	  (*SumsData_)[i].setName(tmpstr);
//	}
//	break;
//      default:
//	std::cerr << "No such variable: " << str << std::endl;
//      }
//    }
//  }

}



void Level1Profiles::adjustPixelIntensities(void)
{
  int len = PixelTimeStamp_.size();
  double tmp = 0.0;
  double maxint = 0.0;
  double factor;

  EventVar<double> avepix;
  valarray<bool> good  = ((PixelTimeStamp_ >= 60.0) &&
			  (PixelTimeStamp_ <= 80.0));

  for(int i=0; i<AvePixVect_.size(); i++)
  {
    avepix = AvePixVect_[i][good];
    tmp    = avepix.max();
    maxint = (maxint > tmp) ? maxint : tmp;
  }

  factor = 1023.0 / maxint;

  for(int i=0; i<AvePixVect_.size(); i++)
  {
    for(int j=0; j<PixelTimeStamp_.size(); j++)
    {
      AvePixVect_[i][j] = AvePixVect_[i][j] * factor;
    }
  }
}



void Level1Profiles::writeExoAtmosphericDataXYZ(void)
{
  int nelem;
  int top[]    = {0,1,2,3,4,5,6};
  int bottom[] = {7,8,9,10,11,12,13};
  int center[] = {24,25,26,27,28,29,30};
  int left[]   = {14,15,16,17,18};
  int right[]  = {19,20,21,22,23};
  int len = PixelTimeStamp_.size();
  int mid = len / 2;
  int first = mid - 200;
  int last  = mid + 200;
  int l;
  double dloc;

  double tmp = 0.0;
  double maxint = 0.0;
  EventVar<double> avepix;
  valarray<bool> good  = ((PixelTimeStamp_ >= 60.0) &&
			  (PixelTimeStamp_ <= 80.0));

  for(int i=0; i<AvePixVect_.size(); i++)
  {
    avepix = AvePixVect_[i][good];
    tmp    = avepix.max();
    maxint = (maxint > tmp) ? maxint : tmp;
  }

//  EventVar<double> ev(Track_Low_X_.size());
//  EventVarVect<double> evv;

//  for(int i=0; i<ev.size(); i++)
//    ev[i] = LocsVect_[3][i] - Track_Low_X_[i];
//  ev.setName("Low X");
//  evv.addEventVar(ev);
//  for(int i=0; i<ev.size(); i++)
//    ev[i] = LocsVect_[10][i] - Track_High_X_[i];
//  ev.setName("High X");
//  evv.addEventVar(ev);
//  evv.Plot("ROI - Track","FPACorrection::Level1Profiles::writeExoAtmosphericDataXYZ()");
  double edge_max_int;
  double extent;
  std::string fname;
  std::string base = GATS_Utilities::ConvertToString(eventNumber_) + ".dat";
  std::ofstream outfile;

  fname = "EdgeData/TopEdge_" + base;
  outfile.open(fname.c_str(),std::ofstream::out);
  outfile.setf(std::ios::scientific, std::ios::floatfield);
  outfile.precision(10);
  edge_max_int =  maxint * XROI_INT / CSUM_INT; 
  nelem = sizeof(top) / sizeof(int);
  for(int i=first; i<=last; i++)
  {
    extent = (Track_High_X_[i] - Track_Low_X_[i]);
    for(int j=0; j<nelem; j++)
    {
      l = top[j];
      dloc = (LocsVect_[l][i] - Track_Low_X_[i]) / extent;
      outfile << dloc << ", ";
      outfile << (AvePixVect_[l][i] / maxint) << "\n";
//      outfile << (AvePixVect_[l][i]) << "\n";
    }
  }
  outfile.close();

  fname = "EdgeData/BottomEdge_" + base;
  outfile.open(fname.c_str(),std::ofstream::out);
  outfile.setf(std::ios::scientific, std::ios::floatfield);
  outfile.precision(10);
  edge_max_int =  maxint * XROI_INT / CSUM_INT; 
  nelem = sizeof(bottom) / sizeof(int);
  for(int i=first; i<=last; i++)
  {
    extent = (Track_High_X_[i] - Track_Low_X_[i]);
    for(int j=0; j<nelem; j++)
    {
      l = bottom[j];
      dloc = (LocsVect_[l][i] - Track_Low_X_[i]) / extent;
      outfile << dloc << ", ";
      outfile << (AvePixVect_[l][i] / maxint) << "\n";
//      outfile << (AvePixVect_[l][i]) << "\n";
    }
  }
  outfile.close();

  fname = "EdgeData/LeftEdge_" + base;
  outfile.open(fname.c_str(),std::ofstream::out);
  outfile.setf(std::ios::scientific, std::ios::floatfield);
  outfile.precision(10);
  edge_max_int =  maxint * YROI_INT / CSUM_INT; 
  nelem = sizeof(left) / sizeof(int);
  for(int i=first; i<=last; i++)
  {
    for(int j=0; j<nelem; j++)
    {
      l = left[j];
      dloc = LocsVect_[l][i] - Track_Low_Y_[i];
      outfile << dloc << ", ";
      outfile << (AvePixVect_[l][i] / maxint) << "\n";
//      outfile << (AvePixVect_[l][i]) << "\n";
    }
  }
  outfile.close();

  fname = "EdgeData/RightEdge_" + base;
  outfile.open(fname.c_str(),std::ofstream::out);
  outfile.setf(std::ios::scientific, std::ios::floatfield);
  outfile.precision(10);
  edge_max_int =  maxint * YROI_INT / CSUM_INT; 
  nelem = sizeof(right) / sizeof(int);
  for(int i=first; i<=last; i++)
  {
    for(int j=0; j<nelem; j++)
    {
      l = right[j];
      dloc = LocsVect_[l][i] - Track_Low_Y_[i];
      outfile << dloc << ", ";
      outfile << (AvePixVect_[l][i] / maxint) << "\n";
//      outfile << (AvePixVect_[l][i]) << "\n";
    }
  }
  outfile.close();

  fname = "EdgeData/Center_" + base;
  outfile.open(fname.c_str(),std::ofstream::out);
  outfile.setf(std::ios::scientific, std::ios::floatfield);
  outfile.precision(10);
  nelem = sizeof(center) / sizeof(int);
  for(int i=first; i<=last; i++)
  {
    extent = (Track_High_X_[i] - Track_Low_X_[i]);
    for(int j=0; j<nelem; j++)
    {
      l = center[j];
      dloc = (LocsVect_[l][i] - Track_Low_X_[i]) / extent;
      outfile << dloc << ", ";
      outfile << (AvePixVect_[l][i] / maxint) << "\n";
//      outfile << (AvePixVect_[l][i]) << "\n";
    }
  }
  outfile.close();
}


void Level1Profiles::dumpData(std::string base, int grid)
{
  std::ofstream outfile;
  std::string fname = base + GATS_Utilities::ConvertToString(eventNumber_) + ".dat";
  outfile.open(fname.c_str(),std::ofstream::out);
  EventVar<double> timestamp;

//  Track_Low_X_.Plot("Track_Low_X_","In dumpData");
//  SolarHighEl_.Plot("SolarHighEl_","In dumpData");

  timestamp = (grid > 0) ? SunSensorTimeStamp_ : PixelTimeStamp_;
  int len = timestamp.size();

  std::cerr << "Track_Low_X_.size() = " << Track_Low_X_.size() << "\n";
  std::cerr << "timestamp.size()    = " << timestamp.size() << "\n";

  for(int i=0;i<len; i++)
  {
    outfile << i << ", " << timestamp[i] << ", ";
    outfile << Track_Low_X_[i] << ", " << Track_High_X_[i] << ", ";
    outfile << Track_Low_Y_[i] << ", " << Track_High_Y_[i] << ", ";
    outfile << ((Track_High_Y_[i]+Track_Low_Y_[i])/2.0) << ", ";
    outfile << ((Track_High_Y_[i]-Track_Low_Y_[i])) << ", ";
    outfile << SolarLowEl_[i] << ", " << SolarHighEl_[i] << ", ";
    outfile << SolarExtent_[i] << ", " << NormalizedSolarExtent_[i] << ", ";
    outfile << DetectorEl_[i]  << ", " << DetectorAz_[i];
    for(int j = 0; j < 31; j++)
      outfile << ", " << AvePixVect_[j][i];
//      outfile << ", " << CenterAvesVect_[j][i];
    outfile << "\n";
  }
  outfile.close();
}



void Level1Profiles::selectSumsData(std::string base, double t1, double t2)
{
  std::ofstream outfile;
  std::string fname = base + GATS_Utilities::ConvertToString(eventNumber_) + ".dat";
  outfile.open(fname.c_str(),std::ofstream::out);
  EventVar<double> timestamp;
  int mafpoints = 25;
  
  int len;
  if(SunSensorTimeStamp_.size() == 0)
    timestamp = PixelTimeStamp_;
  else
    timestamp = SunSensorTimeStamp_;
  len = timestamp.size();
  
  std::valarray<bool> ok = (timestamp > t1) && (timestamp < t2);
  
  for(int i=0;i<len; i++)
  {
    if(ok[i])
    {
      outfile << timestamp[i] << ", " ;
      outfile << CenterAvesVect_[0][i]<< ", ";
      outfile << CenterAvesVect_[1][i]<< ", ";
      outfile << CenterAvesVect_[2][i]<< ", ";
      outfile << CenterAvesVect_[3][i]<< ", ";
      outfile << CenterAvesVect_[4][i]<< ", ";
      outfile << CenterAvesVect_[5][i]<< ", ";
      outfile << CenterAvesVect_[6][i]<< "\n";
    }
  }
  outfile.close();
}









//void Level1Profiles::selectSumsData(std::string base, double t1, double t2)
//{
//  std::ofstream outfile;
//  std::string fname = base + GATS_Utilities::ConvertToString(eventNumber_) + ".dat";
//  outfile.open(fname.c_str(),std::ofstream::out);
//  EventVar<double> timestamp;
//  int mafpoints = 25;
//
//  int len;
//  if(SunSensorTimeStamp_.size() == 0)
//    timestamp = PixelTimeStamp_;
//  else
//    timestamp = SunSensorTimeStamp_;
//  len = timestamp.size();
//  EventVar<double> Ave(len),SmoothAve(len),SD_a(len), SD_sa(len);
//  EventVar<double> s2(len),s3(len),s4(len),s5(len);
//  EventVar<double> sd2(len),sd3(len),sd4(len),sd5(len),sdave(len);
//  EventVar<double> s_sd2(len),s_sd3(len),s_sd4(len),s_sd5(len);
//  EventVar<double> tmp1,tmp2,tmp3,tmp4;
//
//  std::valarray<bool> ok = (timestamp > t1) && (timestamp < t2);
//  std::valarray<bool> ok1 = (timestamp > 30.0) && (timestamp <= 40.0);
//  std::valarray<bool> ok2 = (timestamp > 40.0) && (timestamp <= 50.0);
//  std::valarray<bool> ok3 = (timestamp > 50.0) && (timestamp <= 60.0);
//  std::valarray<bool> ok4 = (timestamp > 60.0) && (timestamp <= 70.0);
//
//  tmp1 = CenterAvesVect_[2];
//  tmp2 = CenterAvesVect_[3];
//  tmp3 = CenterAvesVect_[4];
//  tmp4 = CenterAvesVect_[5];
//
//  Ave  = (tmp1 + tmp2 + tmp3 + tmp4) / 4.0;
//  SmoothAve = Ave.Smooth(mafpoints);
//  SD_a  = Ave.Stddev(mafpoints);
//  SD_sa = SmoothAve.Stddev(mafpoints);
//
//
//  tmp1 = tmp1[ok1];
//  tmp2 = tmp2[ok2];
//  tmp3 = tmp3[ok3];
//  tmp4 = tmp4[ok4];
//
//  len = tmp1.size();
//
//  EventVar<double> ts(len),ave(len),save(len),ave_sd(len),save_sd(len);
//
//  ts = timestamp[ok1];
//  
//  ave = (tmp1 + tmp2 + tmp3 + tmp4) / 4.0;
//  save = ave.Smooth(mafpoints);
//  ave_sd = ave.Stddev(mafpoints);
//  save_sd = save.Stddev(mafpoints);
//
//
//  s2 = CenterAvesVect_[2].Smooth(mafpoints);
//  s3 = CenterAvesVect_[3].Smooth(mafpoints);
//  s4 = CenterAvesVect_[4].Smooth(mafpoints);
//  s5 = CenterAvesVect_[5].Smooth(mafpoints);
//
//  sd2 = CenterAvesVect_[2].Stddev(mafpoints);
//  sd3 = CenterAvesVect_[3].Stddev(mafpoints);
//  sd4 = CenterAvesVect_[4].Stddev(mafpoints);
//  sd5 = CenterAvesVect_[5].Stddev(mafpoints);
//
//  sdave = sd2 + sd3 + sd4 + sd5;
//  sdave = sdave / 4.0;
//
//  s_sd2 = s2.Stddev(mafpoints);
//  s_sd3 = s3.Stddev(mafpoints);
//  s_sd4 = s4.Stddev(mafpoints);
//  s_sd5 = s5.Stddev(mafpoints);
//
//  Ave       = Ave[ok1];
//  SmoothAve = SmoothAve[ok1];
//  SD_a      = SD_a[ok1];
//  SD_sa     = SD_sa[ok1];
//  sdave     = sdave[ok1];
//
//  for(int i=0;i<ts.size(); i++)
//  {
//    outfile << ts[i] << ", " ;
//    outfile << Ave[i] << ", " << SD_a[i] << ", "; 
//    outfile << SmoothAve[i] << ", " << SD_sa[i] << ", ";
//    outfile << ave[i] << ", " << ave_sd[i] << ", "; 
//    outfile << save[i] << ", " << save_sd[i] << ", "; 
//    outfile << sdave[i] << "\n ";
//  }
//  outfile.close();
//
//
////  for(int i=0;i<len; i++)
////  {
////    if(ok[i])
////    {
////      outfile << timestamp[i] << ", " ;
////      outfile << Ave[i] << ", " << SD_a[i] << ", "; 
////      outfile << SmoothAve[i] << ", " << SD_sa[i] << ", ";
////      outfile << CenterAvesVect_[2][i]<< ", " << sd2[i] << ", ";
////      outfile << CenterAvesVect_[3][i]<< ", " << sd3[i] << ", ";
////      outfile << CenterAvesVect_[4][i]<< ", " << sd4[i] << ", ";
////      outfile << CenterAvesVect_[5][i]<< ", " << sd5[i] << ", ";
////      outfile << s2[i]<< ", " << s_sd2[i] << ", ";
////      outfile << s3[i]<< ", " << s_sd3[i] << ", ";
////      outfile << s4[i]<< ", " << s_sd4[i] << ", ";
////      outfile << s5[i]<< ", " << s_sd5[i] << ", ";
////      outfile << sdave[i] << "\n ";
////    }
////  }
////  outfile.close();
//}
//