#include "Level1Data.h"
#include "GATS_DB.h"
#include "CMarkup.h"
#include "GATS_Utilities.hpp"
#include "fortranfunctions.h"


#include <iostream>
#include <cassert>
#include <cmath>

using namespace gatsDBpp;
using namespace GATS_Utilities;

bool Level1Data::AverageList(const std::vector<int>& events, char* l1_ver) {

	std::vector< Vector2D > Signals, SSTran;
	Vector2D TPAlt, Lockdown, SolarExtent, Vangle, RefAngle;
	std::vector<float>  P,T, O_MSIS, O2_MSIS, N2_MSIS;
	std::vector<double>  Exo,Attenuator,ChiSq; //,Lockdown, SolarExtent,Vangle;

    std::vector<float> LLAT, LLON;
    std::vector<int> EEVENT, OORBIT;	
    std::vector<std::string> SSTIME;
	float Curvature=0.;
	float EarthSunD = 0.;

	for(std::vector<int>::const_iterator ievent = events.begin(); ievent != events.end(); ++ievent) {
		if(! Read( *ievent, l1_ver) ) {
			std::cout << " NO data found for event " << *ievent << std::endl;
		} else {

        EEVENT.push_back(*ievent);
        LLAT.push_back(_Lat );
        LLON.push_back(_Lon) ;
        OORBIT.push_back( _OrbitNo);
        SSTIME.push_back( _StartTime);


		Signals.push_back(_DetectorData) ;    // these will get averaged
//		Scans.push_back(_ScanData);
		TPAlt.push_back(_TPAlt) ;                  
		Lockdown.push_back(_LockDownElev);     //these will get averaged
		SolarExtent.push_back(_SolarExtent);
		Vangle.push_back(_ViewAngle);
        RefAngle.push_back(_SSRefraction);
		SSTran.push_back(_SSTransmission);
		
		
        Curvature += _CurvatureRadius;
		EarthSunD += _EarthSunDistance;

/*  altitude is always the same
		if(Z_MSIS.empty() ) Z_MSIS = _msisAlt;
		else std::transform( Z_MSIS.begin(),Z_MSIS.end(), _msisAlt.begin(), Z_MSIS.begin(), std::plus<float>() );
*/	
		if(O_MSIS.empty() ) O_MSIS = _msisO;
		else std::transform( O_MSIS.begin(),O_MSIS.end(), _msisO.begin(), O_MSIS.begin(), std::plus<float>() );

		if(O2_MSIS.empty() ) O2_MSIS = _msisO2;
		else std::transform( O2_MSIS.begin(),O2_MSIS.end(), _msisO2.begin(), O2_MSIS.begin(), std::plus<float>() );

		if(N2_MSIS.empty() ) N2_MSIS = _msisN2;
		else std::transform( N2_MSIS.begin(),N2_MSIS.end(), _msisN2.begin(), N2_MSIS.begin(), std::plus<float>() );
		


/*
		if(Lockdown.empty() ) Lockdown=_LockDownElev;
		else std::transform( Lockdown.begin(), Lockdown.end(), _LockDownElev.begin(), 
			Lockdown.begin(), std::plus<double>() );
	
		if(SolarExtent.empty() ) SolarExtent=_SolarExtent;
		else std::transform( SolarExtent.begin(), SolarExtent.end(), _SolarExtent.begin(), 
			SolarExtent.begin(), std::plus<double>() );
	
		if(Vangle.empty() ) Vangle=_ViewAngle;
		else std::transform( Vangle.begin(), Vangle.end(), _ViewAngle.begin(), 
			Vangle.begin(), std::plus<double>() );
*/
		std::transform(_Pressure.begin(), _Pressure.end(), _Pressure.begin(), log );

		if(P.empty() ) P = _Pressure;
		else std::transform( P.begin(),P.end(), _Pressure.begin(), P.begin(), std::plus<float>() );
	
		if(T.empty() ) T=_Temperature;
		else std::transform( T.begin(),T.end(), _Temperature.begin(), T.begin(), std::plus<float>() );

/*
for(unsigned h=0; h< _Pressure.size(); ++h) 

std::cout << _Temperature[h] << "  " << _Pressure[h] << std::endl;

std::cout <<"&" << std::endl;
*/
	
		if(Exo.empty() ) Exo=_exoV;
		else std::transform( Exo.begin(),Exo.end(), _exoV.begin(), Exo.begin(), std::plus<double>() );
	
		if(Attenuator.empty() ) Attenuator=_Attenuators;
		else std::transform( Attenuator.begin(), Attenuator.end(), _Attenuators.begin(), 
			Attenuator.begin(), std::plus<double>() );
// NO Chi Sq
        if(ChiSq.empty() ) ChiSq=_ThermalChiSq;
        else std::transform( ChiSq.begin(), ChiSq.end(), _ThermalChiSq.begin(),
            ChiSq.begin(), std::plus<double>() );
	
		}

	}
	
	if(Signals.empty() ) {
		std::cout  << "NO signals found to average " << std::endl ;
		return false;
	}

	int num_events = (int)Signals.size();
	
	std::transform( P.begin(),P.end(), P.begin(), std::bind2nd( std::divides<float>(),(float)num_events ) );
	std::transform( T.begin(),T.end(), T.begin(), std::bind2nd( std::divides<float>(),(float)num_events ) );
/*
for(unsigned h=0; h< P.size(); ++h) 

std::cout << T[h] << " " << P[h] << std::endl;

std::cout <<"&" << std::endl;
exit(23);
*/
	std::transform( O_MSIS.begin(),O_MSIS.end(), O_MSIS.begin(),std::bind2nd( std::divides<float>(),(float)num_events )  );
	std::transform( O2_MSIS.begin(),O2_MSIS.end(), O2_MSIS.begin(),std::bind2nd( std::divides<float>(),(float)num_events )  );	
	std::transform( N2_MSIS.begin(),N2_MSIS.end(), N2_MSIS.begin(),std::bind2nd( std::divides<float>(),(float)num_events )  );
	
	std::transform( Exo.begin(),Exo.end(), Exo.begin(), std::bind2nd( std::divides<double>(),(double)num_events ) );
	std::transform( Attenuator.begin(),Attenuator.end(), Attenuator.begin(), 
		std::bind2nd( std::divides<double>(),(double)num_events ) );
    std::transform( ChiSq.begin(),ChiSq.end(), ChiSq.begin(),
        std::bind2nd( std::divides<double>(),(double)num_events ) );

// start signal averaging
// first we need the data range that covers all the events
	double zmin = -10. , zmax = 900.;
	std::vector<double> g;
	std::vector<double>::iterator iter;

	for(int i = 0; i < num_events; ++i ) {
		g = TPAlt[i];		
	 	iter = std::min_element(g.begin(), g.end() );
		zmin = std::max( *iter, zmin) ;
		iter = std::max_element(g.begin(), g.end() );
		zmax = std::min( *iter, zmax);
//std::cout << " event " << i << " zmin " << zmin << "  zmax " << zmax << "  " << g.size() << std::endl;
	}
	
//	exit(23);
	

	
	unsigned u1, u2;

	for(int i = 0; i < num_events; ++i ) {
		g = TPAlt[i];

		iter = std::find(g.begin(), g.end(), zmax) ;
		u1= std::distance(g.begin(), iter);
		iter = std::find(g.begin(), g.end(), zmin) ;
		u2 = std::distance(g.begin(), iter)+1;

		if( i == 0 ) { 
			_TPAlt = std::vector<double>(g.begin()+u1, g.begin()+u2);
		}

		Vector2D Sigs = Signals[i];

		for(unsigned j = 0; j < Sigs.size(); ++j) {
			g = Sigs[j];
			g = std::vector<double>(g.begin()+u1, g.begin()+u2);
			if( i == 0 ) _DetectorData[j] = g;
			else std::transform(_DetectorData[j].begin(), _DetectorData[j].end(), g.begin(), 
				_DetectorData[j].begin(), std::plus<double>() );
				

/*
if(j==2) {
for(unsigned h=0; h< g.size(); ++h) 

std::cout << g[h] << "  " << _TPAlt[h] << std::endl;

std::cout <<"&" << std::endl;
exit(23);
}
*/
		}
		
// lockdown
		g = Lockdown[i];
		g = std::vector<double>(g.begin()+u1, g.begin()+u2);
		
		if( i == 0 ) _LockDownElev = g;
		else std::transform(_LockDownElev.begin(), _LockDownElev.end(), g.begin(), 
			_LockDownElev.begin(), std::plus<double>() );		

		
// solar extent

		g = SolarExtent[i];
		g = std::vector<double>(g.begin()+u1, g.begin()+u2);
		
		if( i == 0 ) _SolarExtent = g;
		else std::transform(_SolarExtent.begin(), _SolarExtent.end(), g.begin(), 
			_SolarExtent.begin(), std::plus<double>() );
// vangle

		g = Vangle[i];
		g = std::vector<double>(g.begin()+u1, g.begin()+u2);
		
		if( i == 0 ) _ViewAngle = g;
		else std::transform(_ViewAngle.begin(), _ViewAngle.end(), g.begin(), 
			_ViewAngle.begin(), std::plus<double>() );

// REfraction
        g = RefAngle[i];
        g = std::vector<double>(g.begin()+u1, g.begin()+u2);

 	   	if( i == 0 ) _SSRefraction = g;
    	else std::transform(_SSRefraction.begin(), _SSRefraction.end(), g.begin(),
            _SSRefraction.begin(), std::plus<double>() );

// sunsensor transmission
//   SStran now contains the composite average transmission and the individual transmissions from ROI's 1-7


        Sigs = SSTran[i];

        for(unsigned j = 0; j < Sigs.size(); ++j) {
            g = Sigs[j];
            g = std::vector<double>(g.begin()+u1, g.begin()+u2);
  Sigs[j] =g;
            if( i == 0 ) _SSTransmission[j] = g;
            else std::transform(_SSTransmission[j].begin(), 
                                _SSTransmission[j].end(), g.begin(),
                                _SSTransmission[j].begin(), std::plus<double>() );
        }
 SSTran[i] = Sigs;



    } // end loop over events

// variance of sunsensor transmission

      for(int i = 0; i<num_events; ++i) {

        Vector2D Sigs = SSTran[i];
	
	


     /*  
Vector2D oneAlt;
        for(int j = 0; j< Sigs.size() ; ++j)   {

 oneAlt[j].push_back(S


        }
*/

     }
//exit(23);

//   now we need to the mean and standard dev of the transmittance


/*

double sum = std::accumulate(std::begin(v), std::end(v), 0.0);
double m =  sum / v.size();

double accum = 0.0;
std::for_each (std::begin(v), std::end(v), [&](const double d) {
    accum += (d - m) * (d - m);
});

double stdev = sqrt(accum / (v.size()-1));    	

*/


	for(unsigned j=0; j< _DetectorData.size(); ++j) {
		std::transform( _DetectorData[j].begin(), _DetectorData[j].end(), 
			_DetectorData[j].begin(), 
			std::bind2nd( std::divides<double>(),(double)num_events ) );

      std::transform( _DetectorData[j].begin(), _DetectorData[j].end(),
         _DetData[j].begin(), std::bind2nd( std::multiplies<double>(),29000.00 ) );

/*
if(j==1) {
for(unsigned h=0; h< _TPAlt.size(); ++h) 

std::cout << _DetectorData[j][h] << "  " << _TPAlt[h] << std::endl;

std::cout <<"&" << std::endl;
exit(23);
}
*/
	}
	
	std::transform( _LockDownElev.begin(), _LockDownElev.end(), 
			_LockDownElev.begin(), 
			std::bind2nd( std::divides<double>(),(double)num_events ) );

	std::transform( _SolarExtent.begin(), _SolarExtent.end(), 
			_SolarExtent.begin(), 
			std::bind2nd( std::divides<double>(),(double)num_events ) );
			
	std::transform( _ViewAngle.begin(), _ViewAngle.end(), 
			_ViewAngle.begin(), 
			std::bind2nd( std::divides<double>(),(double)num_events ) );


    std::transform( _SSRefraction.begin(), _SSRefraction.end(),
            _SSRefraction.begin(),
            std::bind2nd( std::divides<double>(),(double)num_events ) );


    for(unsigned j=0; j< _SSTransmission.size(); ++j) 
        std::transform( _SSTransmission[j].begin(), _SSTransmission[j].end(),
            _SSTransmission[j].begin(),
            std::bind2nd( std::divides<double>(),(double)num_events ) );

/*
for(unsigned h=0; h< _TPAlt.size(); ++h) 

std::cout << ConvertToStringPrec(_SSRefraction[h]) << "  " << _TPAlt[h] << std::endl;

std::cout <<"&" << std::endl;
exit(23);
*/


// end of signal averaging


std::transform(P.begin(), P.end(), P.begin(), exp );
	_Pressure = P;
	_Temperature = T;

	_exoV = std::vector<double>(_exoV.size(), 29000.00) ; //Exo;
	_Attenuators = Attenuator;
    _ThermalChiSq = ChiSq;
	_msisO = O_MSIS;
	_msisO2 = O2_MSIS;
	_msisN2 = N2_MSIS;
	_CurvatureRadius = Curvature / num_events;
 	_EarthSunDistance = EarthSunD / num_events;

// pick middle val or near the middle
 if (  num_events % 2 ) num_events++;
 num_events /= 2;
 num_events -= 1;
      _EventNo=   EEVENT.at(num_events);
    _Lat  =  LLAT.at(  num_events );
    _Lon  =  LLON.at(num_events); 
     _OrbitNo =  OORBIT.at( num_events);
    _StartTime =   SSTIME.at( num_events);



//std::cout << "done " << std::endl;
//exit(23);
	return true;
}



bool Level1Data::Read(const int EventNo, char* l1_ver)
{

	if( ! ReadSummary( EventNo,l1_ver) ) { 
		std::cout << " Level 1 Data not found for Event " << EventNo 
			<< " Version: " << l1_ver << std::endl;
		return false; 
	}  // Event does not exist

	ReadScans(EventNo);
	ReadZPT(EventNo);
	ReadGridded(EventNo);

	ReadDriftData(EventNo);		
	ReadEphem(EventNo);
	ReadMSIS(EventNo);
	
	return true;	
}

void Level1Data::Init(CMarkup* xml,const std::string& section)
{
	xml->InsideRoot();
	if(xml->FindElem( section.c_str() ) ) {
		xml->IntoElem() ;
		assert(xml->FindElem("dbConnection") );
		std::string value = CMarkup::TrimString( xml->GetData() ); 
		DatabaseHandles* handles = GetDatabaseHandles();
		_dbconn = handles->Get(value);
		assert(_dbconn);
		
		xml->ResetMainPos();
		xml->ResetMainPos();		
		_summaryTable = (xml->FindElem("summaryTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "l1_sofie_summary" ;
			 
		xml->ResetMainPos();	 
		_mergedTPtable = (xml->FindElem("mergeTPTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "MERGED_PROFILES" ;
			 
		xml->ResetMainPos();	 
		_transmissionTable = (xml->FindElem("transmissionTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_TRANSMISSIONS" ;

		xml->ResetMainPos();
		_signalTable = (xml->FindElem("signalTable") ) ?
   			CMarkup::TrimString( xml->GetData() ) : "GRIDDED_SIGNAL" ;


		xml->ResetMainPos();	 
		_statsTable = (xml->FindElem("statsTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "DETECTOR_STATS" ;
			 
		xml->ResetMainPos();	 
		_solarAngleTable = (xml->FindElem("solarAngleTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "SOLAR_SOURCE_ANGLE" ;
			 
		xml->ResetMainPos();	 
		_solarTable = (xml->FindElem("solarTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "SOLAR_SOURCE" ;
		
		xml->ResetMainPos();	 
		_ephemTable = (xml->FindElem("ephemerisTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_EPHEMERIS" ;
			 	 
		xml->ResetMainPos();	 
		_timeTable = (xml->FindElem("timeTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_TIME" ;

		xml->ResetMainPos();	 
		_lockdownTable = (xml->FindElem("lockdownTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_SOLAR_VIEW" ;

		xml->ResetMainPos();	 
		_altsTable = (xml->FindElem("altsTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_IND_VAR" ;
			 
		xml->ResetMainPos();	 
		_driftTable = (xml->FindElem("driftTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "DRIFT_CORRECTIONS" ;

		xml->ResetMainPos();	 
		_solarDatTable = (xml->FindElem("solarViewTable") ) ? 
			 CMarkup::TrimString( xml->GetData() ) : "GRIDDED_SOLAR_VIEW" ;

        xml->ResetMainPos();
        _miscTable = (xml->FindElem("miscellaneousTable") ) ?
             CMarkup::TrimString( xml->GetData() ) : "MISCELLANEOUS" ;

        xml->ResetMainPos();
        _msisTable = (xml->FindElem("MSISTable") ) ?
		     CMarkup::TrimString( xml->GetData() ) : "MSIS_XRATIOS" ;

		xml->OutOfElem();
	}
}

bool Level1Data::ReadSummary(const int EventNo, char* l1_ver)
{
	std::string sql = std::string("select ") 
		+ "Event,orbit_number,mode,latitude,longitude, start_time,end_time "
		+ " from " + _summaryTable 
		+ " where Event = " + ConvertToString(EventNo)
                + " and L1_Version = '" + l1_ver + "'";
	
	GATS_DB_Results R = _dbconn->Query(sql.c_str() );
	if( R.empty() ) { return false ; }
	
	_EventNo = R[0][0].asUlong();
	_OrbitNo = R[0][1].asUlong();
	_Mode = R[0][2].asString();
	_Lat = R[0][3].asFloat();
	_Lon = R[0][4].asFloat();
	_StartTime = R[0][5].asString();
	_EndTime =   R[0][6].asString();
	
	return true;
}

void Level1Data::ReadScans(const int  EventNo)
{
	GATS_DB_Results R;
	std::string sql;
	
	_ElevAng.clear();
	_ScanData.clear();
	
	sql = std::string("select Array") 
		+" from " + _solarAngleTable 
		+ " where Event = " + ConvertToString(EventNo) 
		+ " AND Parameter = 'Angle'";

	R = _dbconn->Query(sql.c_str() );
	assert(!R.empty() );
	_ElevAng = R[0][0].asVecDouble();
				
	sql = "select Array from "  + _solarTable 
		+ " where Event = " + ConvertToString(EventNo) 
		+ " order by Signal_no";
			
	R = _dbconn->Query(sql.c_str() );
	assert((int)R.size() == 16);
			
	for (GATS_DB_Results::const_iterator resIter = R.begin() ; 
		resIter != R.end(); ++resIter) {
			_ScanData.push_back( (*resIter)[0].asVecDouble() );
	}	
// Add a dummy for the sunsensor this would be for the sunsensor ozone

    _ScanData.push_back(std::vector<double>(_ElevAng.size(), 1.0) );

// sanity checks
	assert(! _ElevAng.empty() );
	assert(monotonic_if( _ElevAng.begin(), _ElevAng.end(), std::less<double>()    )  );
	for(unsigned i = 0 ; i < _ScanData.size();  ++i) {
		assert(_ScanData[i].size() == _ElevAng.size() );	
	}
	
}

void Level1Data::ReadGridded(const int  EventNo)
{
	GATS_DB_Results R;
	std::string sql;

	_DetectorData.clear();
    _DetData.clear();
	_Attenuators.clear();
    _ThermalChiSq.clear();
	_FilterTemps.clear();
	_Times.clear();
	_TPAlt.clear();
	_ViewAngle.clear();
	_LockDownElev.clear();
	_SolarExtent.clear();
	_SSRefraction.clear();
    _SSTransmission.clear();
	

	sql = "select Array from " + _altsTable 
			+ " where Event = " + ConvertToString(EventNo) 
			+ " And Parameter = 'Altitude'" ;	
	R = _dbconn->Query(sql.c_str() );
	assert( !R.empty() );
	_TPAlt = R[0][0].asVecDouble();
	std::reverse(_TPAlt.begin(), _TPAlt.end() );
		
	sql = "select Array from " + _altsTable 
			+ " where Event = " + ConvertToString(EventNo) 
			+ " And Parameter = 'Time'" ;
	R = _dbconn->Query(sql.c_str() );
	assert( !R.empty() );
	_Times = R[0][0].asVecDouble();
	std::reverse(_Times.begin(), _Times.end() );
			
	sql = "select Array from " + _altsTable 
			+ " where Event = " + ConvertToString(EventNo) 
			+ " And Parameter = 'Angle'";
	R = _dbconn->Query(sql.c_str() );
	assert( !R.empty() );
	_ViewAngle = R[0][0].asVecDouble();
	std::reverse(_ViewAngle.begin(), _ViewAngle.end() );
		
	sql = "select Array from " + _transmissionTable 
			+ " where Event = " + ConvertToString(EventNo) 
			+ " order by Signal_No" ;
	R = _dbconn->Query(sql.c_str() );
	assert((int)R.size() == 24);
	for (GATS_DB_Results::const_iterator resIter = R.begin(); 
		resIter != R.end(); ++resIter) {
			_DetectorData.push_back( (*resIter)[0 ].asVecDouble() );
		std::reverse( (_DetectorData.back()).begin(), (_DetectorData.back()).end() );
	}

    sql = "select Array from " + _signalTable
            + " where Event = " + ConvertToString(EventNo)
            + " order by Signal_No" ;
    R = _dbconn->Query(sql.c_str() );
    assert((int)R.size() == 24);

    for (GATS_DB_Results::const_iterator resIter = R.begin();
        resIter != R.end(); ++resIter) {
            _DetData.push_back( (*resIter)[0 ].asVecDouble() );
		std::reverse( (_DetData.back()).begin(), (_DetData.back()).end() );
     }

	sql = "select Array from " + _solarDatTable
		+ " where Event = " + ConvertToString(EventNo) 
	  	+ " And Parameter = 'EffLockdownEl'"; //'LockdownEl'" ;
	R = _dbconn->Query(sql.c_str() );
	assert(!R.empty() );
	_LockDownElev = R[0][0].asVecDouble() ;
	std::reverse(_LockDownElev.begin(), _LockDownElev.end() );

	sql = "select Array from " + _solarDatTable
  		+ " where Event = " + ConvertToString(EventNo) 
  		+ " And Parameter = 'SolarExtent'" ;
	R = _dbconn->Query(sql.c_str() );
	assert(!R.empty() );
	_SolarExtent= R[0][0].asVecDouble() ;
	std::reverse(_SolarExtent.begin(), _SolarExtent.end() );


	sql = "select Array from " + _solarDatTable
  		+ " where Event = " + ConvertToString(EventNo) 
  		+ " And Parameter = 'RefractionAngle'" ;
	R = _dbconn->Query(sql.c_str() );
	if(!R.empty() )  {
		_SSRefraction = R[0][0].asVecDouble() ;
		std::reverse(_SSRefraction.begin(), _SSRefraction.end() );
	}
// sunsensor transmission


    sql = "select Array from " + _solarDatTable
       + " where Event = " + ConvertToString(EventNo)
       + " And Parameter = 'SSTransmission'" ;
    R = _dbconn->Query(sql.c_str() );
    if(!R.empty() )  {
       _SSTransmission.push_back( R[0][0].asVecDouble() );
       std::reverse(_SSTransmission[0].begin(), _SSTransmission[0].end() );
     }


/*
     for(int iRegion = 1; iRegion <=7; ++iRegion) {
         sql = "select Array from GRIDDED_AVG_PIXELS where Event = " + ConvertToString(EventNo)
             + " And Signal_no = "+ConvertToString(iRegion)  ;
         R = _dbconn->Query(sql.c_str() );
         assert(!R.empty() ); 
         _SSTransmission.push_back( R[0][0].asVecDouble()) ;
         std::reverse(_SSTransmission[iRegion].begin(), _SSTransmission[iRegion].end() );
     }
*/

// get some other data
		sql =std::string( "select Value  from " )
		+ _statsTable 
		+ " where Event = " + ConvertToString(EventNo) 
		+ " And Parameter = 'Attenuation'"
		+ " Order by Signal_No";
		
		R = _dbconn->Query(sql.c_str() );
		assert( R.size() ==16);
		for (GATS_DB_Results::const_iterator resIter = R.begin(); 
		resIter != R.end(); ++resIter) {
			_Attenuators.push_back( (*resIter)[0].asDouble() );
		}
		
		sql =std::string( "select Value+273.15  from " )
		+ _statsTable 
		+ " where Event = " + ConvertToString(EventNo) 
		+ " And Parameter = 'FilterTemperature'"
		+ " Order by Signal_No";
		
		R = _dbconn->Query(sql.c_str() );
		assert( R.size() ==16);
		for (GATS_DB_Results::const_iterator resIter = R.begin(); 
		resIter != R.end(); ++resIter) {
			_FilterTemps.push_back( (*resIter)[0].asDouble() );
		}
// New NO Thermal correction
        sql =std::string( "select Value  from " )
        + _statsTable
        + " where Event = " + ConvertToString(EventNo)
        + " And Parameter = 'ThermalChiSq'"
        + " Order by Signal_No";

        R = _dbconn->Query(sql.c_str() );
        assert( R.size() ==16);
        for (GATS_DB_Results::const_iterator resIter = R.begin();
        resIter != R.end(); ++resIter) {
            _ThermalChiSq.push_back( (*resIter)[0].asDouble() );
        }


// let's filter out the missing values
	std::vector<double>::iterator res = 
		std::find_if(_TPAlt.begin(), _TPAlt.end(), std::bind2nd(std::greater<double>(),0.) );
	unsigned ik = std::distance(_TPAlt.begin(), res);
	if(ik != 0) {
		_TPAlt.erase( _TPAlt.begin(), _TPAlt.begin()+ik );		
		_Times.erase( _Times.begin(), _Times.begin()+ik );
		_ViewAngle.erase( _ViewAngle.begin(), _ViewAngle.begin()+ik );
		_LockDownElev.erase( _LockDownElev.begin(), _LockDownElev.begin()+ik );
		_SolarExtent.erase( _SolarExtent.begin(), _SolarExtent.begin()+ik );
        _SSRefraction.erase( _SSRefraction.begin(), _SSRefraction.begin()+ik );
//        _SSTransmission.erase( _SSTransmission.begin(), _SSTransmission.begin()+ik );
		for(unsigned i=0; i< _DetectorData.size(); ++i) {
			_DetectorData[i].erase( _DetectorData[i].begin(), _DetectorData[i].begin()+ik   )   ;
		}
        for(unsigned i=0; i< _DetData.size(); ++i) {
            _DetData[i].erase( _DetData[i].begin(), _DetData[i].begin()+ik   )   ;
        }

        for(unsigned i = 0; i <=_SSTransmission.size() ; ++i) {
           _SSTransmission[i].erase( _SSTransmission[i].begin(), _SSTransmission[i].begin()+ik  );
        }  

	}

	res = std::find_if(_TPAlt.begin(), _TPAlt.end(), std::bind2nd(std::less<double>(),0.) );
	ik = std::distance(_TPAlt.begin(), res);
	if(ik != _TPAlt.size() ) {
		_TPAlt.erase( _TPAlt.begin()+ik,  _TPAlt.end() );		
		_Times.erase( _Times.begin()+ik,  _Times.end()  );
		_ViewAngle.erase(  _ViewAngle.begin()+ik,  _ViewAngle.end()  );
		_LockDownElev.erase(  _LockDownElev.begin()+ik, _LockDownElev.end() );
		_SolarExtent.erase(  _SolarExtent.begin()+ik , _SolarExtent.end());
        _SSRefraction.erase(  _SSRefraction.begin()+ik , _SSRefraction.end());
//        _SSTransmission.erase(  _SSTransmission.begin()+ik , _SSTransmission.end());
		for(unsigned i=0; i< _DetectorData.size(); ++i) {
			_DetectorData[i].erase(  _DetectorData[i].begin()+ik,_DetectorData[i].end()   )   ;
		}

        for(unsigned i=0; i< _DetData.size(); ++i) {
            _DetData[i].erase(  _DetData[i].begin()+ik,_DetData[i].end()   )   ;
        }
        for(unsigned i = 0; i <=_SSTransmission.size() ; ++i) {
           _SSTransmission[i].erase( _SSTransmission[i].begin()+ik, _SSTransmission[i].end()  );
        }

	}
// now some sanity checks

	 res =std::find_if(_Attenuators.begin(), _Attenuators.end(), 
		std::bind2nd(std::less<double>(),0.) );
	assert(res == _Attenuators.end() );
	res =std::find_if(_Attenuators.begin(), _Attenuators.end(),
                std::bind2nd(std::greater<double>(),1.) );
	assert(res == _Attenuators.end() );

	assert( ! _TPAlt.empty() );
	assert(monotonic_if( _TPAlt.begin(), _TPAlt.end(), std::greater<double>()    )  );
	assert( _Times.size() == _TPAlt.size() && 
           _ViewAngle.size() == _TPAlt.size() &&
		   _LockDownElev.size() == _TPAlt.size() && 
          _SolarExtent.size() == _TPAlt.size() &&
          _SSRefraction.size() == _TPAlt.size() );
       //   _SSTransmission.size() == _TPAlt.size()  );
	for(unsigned i=0; i< _DetectorData.size(); ++i) {
		assert(_DetectorData[i].size() == _TPAlt.size() );
	}
    for(unsigned i=0; i< _DetData.size(); ++i) {
        assert(_DetData[i].size() == _TPAlt.size() );
    }

   for(unsigned i = 0; i < _SSTransmission.size() ; ++i) {
       assert(_SSTransmission[i].size() == _TPAlt.size() );
   }
	
}

void Level1Data::ReadZPT(const int EventNo)
{
	GATS_DB_Results R;
	std::string sql;

	_Altitude.clear();
	_Pressure.clear();
	_Temperature.clear();

	std::vector<double> fred;
	sql = "select Array from "+ _mergedTPtable +
		+ " where Event = " + ConvertToString(EventNo) 
		+ " And Parameter = 'Altitude'";
	
	R = _dbconn->Query(sql.c_str() );
	assert( !R.empty() );
	fred = R[0][0].asVecDouble();
	_Altitude = std::vector<float>(fred.begin(),fred.end());


	sql = "select Array from "+ _mergedTPtable +
		+ " where Event = " + ConvertToString(EventNo) 
		+ " And Parameter = 'Temperature'";
	
	R = _dbconn->Query(sql.c_str() );
	fred = R[0][0].asVecDouble();
	_Temperature = std::vector<float>(fred.begin(),fred.end());
		
	sql = "select Array from "+ _mergedTPtable +
		+ " where Event = " + ConvertToString(EventNo) 
		+ " And Parameter = 'Pressure'";

	R = _dbconn->Query(sql.c_str() );
	fred = R[0][0].asVecDouble();
	_Pressure = std::vector<float>(fred.begin(),fred.end());

	assert(!_Altitude.empty() && _Pressure.size() == _Altitude.size() 
       && _Temperature.size() == _Altitude.size() );	
	assert(monotonic_if( _Altitude.begin(), _Altitude.end(), std::greater<float>()  )  );
}

void Level1Data::ReadMSIS(const int EventNo)
{
	_msisAlt.clear();
	_msisO.clear();
	_msisO2.clear();
	_msisN2.clear();

	if(! _dbconn->TableNames( _msisTable.c_str() ).empty() ) {
		
		std::string sql = "select Array from "
		+ _msisTable + " where Parameter = 'Altitude' And Event ="+ ConvertToString(EventNo);
		GATS_DB_Results R = _dbconn->Query(sql.c_str() );
		if(!R.empty())  {
			std::vector<double> T = R[0][0].asVecDouble();
			_msisAlt = std::vector<float>(T.size());
			std::reverse_copy( T.begin(), T.end(), _msisAlt.begin() );

			sql = "select Array from "
			   + _msisTable + " where Parameter = 'O_MixRatio' And Event ="+ ConvertToString(EventNo);
			R = _dbconn->Query(sql.c_str() );
			T= R[0][0].asVecDouble();  _msisO = std::vector<float>(T.size());
			std::reverse_copy( T.begin(), T.end(), _msisO.begin() );
//std::cout << "Cut msis O factor of 10 " << std::endl;
//std::transform( _msisO.begin(),_msisO.end(), _msisO.begin(), std::bind2nd( std::divides<float>(),10.0 ) );


			sql = "select Array from "
			   + _msisTable + " where Parameter = 'O2_MixRatio' And Event ="+ ConvertToString(EventNo);
			R = _dbconn->Query(sql.c_str() );
			R = _dbconn->Query(sql.c_str() );
			T= R[0][0].asVecDouble();  _msisO2 = std::vector<float>(T.size());
			std::reverse_copy( T.begin(), T.end(), _msisO2.begin() );


			sql = "select Array from "
			   + _msisTable + " where Parameter = 'N2_MixRatio' And Event ="+ ConvertToString(EventNo);
			R = _dbconn->Query(sql.c_str() );
			R = _dbconn->Query(sql.c_str() );
			T= R[0][0].asVecDouble();  _msisN2 = std::vector<float>(T.size());
			std::reverse_copy( T.begin(), T.end(), _msisN2.begin() );
		}
		assert(_msisAlt.size() == _msisO.size() && 
			_msisAlt.size() == _msisO2.size() &&
			_msisAlt.size() == _msisN2.size() );
/*
for(unsigned j = 0; j < _msisAlt.size(); j++ ) {
		std::cout << _msisAlt[j] << " " << _msisO[j] << " " << _msisO2[j] << " " << _msisN2[j] << std::endl;
}
exit(23);
*/
	}


}

void Level1Data::ReadEphem(const int EventNo)
{
	std::string sql;
	GATS_DB_Results R;
	_TPLat = std::vector<double>(_DetectorData[0].size(), _Lat); 
	_TPLon = std::vector<double>(_DetectorData[0].size(), _Lon); 
// use radgrav to get radius
//static float zero=0.0;
//        float grav;
//	::radgrav_( &zero, &_Lat, &grav, &_CurvatureRadius);

	if(! _dbconn->TableNames( _miscTable.c_str() ).empty() ) {
		
		sql = "select Value from "+ _miscTable + " where Event ="+ ConvertToString(EventNo)
			+ " and Parameter = 'RadiusOfCurvature'";
		R = _dbconn->Query(sql.c_str() );
		if(!R.empty()) {
			 _CurvatureRadius = R[0][0].asFloat();
		}

		sql = "select Value from "+ _miscTable + " where Event ="+ ConvertToString(EventNo)
			+ " and Parameter = 'BearingLOS'";
		R = _dbconn->Query(sql.c_str() );
		if(!R.empty()) {
			 float V =  R[0][0].asFloat();
// Convert from Radians to Degrees
			_LOS_Bearing = ((V < 0.0) ? 2.0* M_PI+V : V) * 180./M_PI ; 
		}

		sql = "select Value from "+ _miscTable + " where Event ="+ ConvertToString(EventNo)
			+ " and Parameter = 'EarthSunDistance'";
		R = _dbconn->Query(sql.c_str() );
		if(!R.empty()) {
			 _EarthSunDistance = R[0][0].asFloat()  ;
		}
		
	}
}


void Level1Data::ReadDriftData(const int EventNo)
{
	_exoV.clear();
	_noise.clear();
	std::string sql =std::string( "select Value  from " )
	  + _driftTable 
	  + " where Event = " + ConvertToString(EventNo) 
	  + " And Parameter = 'SignalRef'"
	  + " Order by Signal_No";
	  
	GATS_DB_Results R = _dbconn->Query(sql.c_str() );
	assert( R.size() >15 );
	for (GATS_DB_Results::const_iterator resIter = R.begin(); 
	resIter != R.end(); ++resIter) {
		_exoV.push_back( (*resIter)[0].asDouble() );
	}

	sql =std::string( "select Value  from " )
	  + _driftTable 
	  + " where Event = " + ConvertToString(EventNo) 
	  + " And Parameter = 'StDev'"
	  + " Order by Signal_No";
	  
	R = _dbconn->Query(sql.c_str() );
	assert( R.size() >15 );
	for (GATS_DB_Results::const_iterator resIter = R.begin(); 
	resIter != R.end(); ++resIter) {
		_noise.push_back( (*resIter)[0].asDouble() );
	}	
}