//
//  $Id$
//-----------------------------------------------------------------------
//
//            (c) Copyright 2006 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.
//
//-----------------------------------------------------------------------
//
// Module:	SolarSourceModel.cpp
//
// Author:	John Burton	<john@gats-inc.com>
//
// Date:  	Thu May  4 17:01:31 2006
//
//-----------------------------------------------------------------------
//
// Modification History:
//
//     $Log$
//
//-----------------------------------------------------------------------
//
//-----------------------------------------------------------------------
// Include Files:
//-----------------------------------------------------------------------
//

#include "ConfigFile.h"
#include "Event.h"
#include "SolarSourceModel.h"
#include "AllenCurve.h"
#include "GATS_Utilities.hpp"

#include <iostream>
#include <cassert>
//
//-----------------------------------------------------------------------
// Defines and Macros:
//-----------------------------------------------------------------------
//



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



//
//-----------------------------------------------------------------------
// Utility Routines:
//-----------------------------------------------------------------------
//
using namespace GATS_Utilities;

int SolarSourceModel(Event& L0, Event& L1, Event& Tmp, Event& SD, ConfigFile& cf)
{
// SOFIE Band 1 through Band 16 wavelengths.

double band_waveNum[] =
        {
        0.291, // modification to the AllenCurve class allows this wavelength
        0.330, 0.867, 1.04, 2.46, 2.62, 2.79, 2.94,
          3.06,  3.12,  3.38,  3.48, 4.32, 4.65, 5.01, 5.32 };

if( cf.ValidEntry("SolarSourceModel","Wavelengths")   ) {
	std::vector<float> waveN = cf.GetRealValues("SolarSourceModel","Wavelengths");
	assert((int)waveN.size() == 16 );
	std::cout <<  "Replace default wavelengths for solar source model" << std::endl;
    std::copy(waveN.begin(), waveN.end(), band_waveNum);    

}	
	EventVarVect<double> Evv("L1_SolarModel", "Normalized 1-D SLDC for all 16 SOFIE Bands");
	
	std::string distanceId("L1_SolarDistance");
	double solarDistance = 150e6;   // 150 million km
	if( L1.EventVarExists(distanceId) ){
		EventVar<double> Ev;
		L1.getEventVar("L1_SolarDistance",Ev);
		solarDistance = Ev[0];
	} else {
		cerr << distanceId << " not found using default solar distance " << std::endl;
	}

//create an arbitrary angle vector for the solar model in Arcminutes;   
// this angle is relative to the solar top edge.  Negative angles are above the Sun,
// positive angles are down or Earthward on the solar disk.
// 
// this should probably be setup from configuration options.
	vector<double> el_sol;
	el_sol.push_back(-3.5);
    while( el_sol.back() < 34.7) el_sol.push_back(el_sol.back()+0.14);

	AllenCurve allenCurve;
	
	for(int iBand = 0; iBand < 16; ++iBand) {
		vector<double> sldc;
		for(vector<double>::const_iterator position=el_sol.begin() ; 
		position != el_sol.end() ; ++position) {		
			sldc.push_back( allenCurve.getSolarIntensity(*position, 
				band_waveNum[iBand], solarDistance ) );				
		}  
		std::string name = string("L1_SolarCurveBand_")+ConvertToString(iBand+1);
		std::string descrip = string("Normalized 1-D SLDC for SOFIE Band ")
			+ ConvertToString(iBand+1) + string(" wavelength ") + ConvertToString(band_waveNum[iBand] )
			+ string(" microns");
//std::cout << descrip << std::endl;	
		EventVar<double> Ev( name, descrip, "Normalized Response", &sldc[0], sldc.size() );
		Ev.setAxisRange(0.f, 1.1f);	
		Evv.addEventVar( Ev );			
	}
//exit(23);
    L1.addEventVar(Evv);
    	
//convert lockdown from armin to radians	
	transform(el_sol.begin(), el_sol.end(), el_sol.begin(),
		std::bind2nd(std::multiplies<double>(), M_PI/10800.f) );  //convert to radians
		
//	EventVar<double> Ev("L1_SLDCAngleFromTopEdge", "Angle from Top Edge of SLDC", "radians", &el_sol[0], el_sol.size() );
//	Evv.addEventVar(Ev);      
//	L1.addEventVar(Evv);
    
	EventVar<double> Ev("L1_SLDCAngleFromTopEdge", "Angle from Top Edge of SLDC", "radians", &el_sol[0], el_sol.size() );
	L1.addEventVar(Ev);
    
	return 0;
}