//
//  $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:	RefractionAngleCalc.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 "EventVar.h"
#include "RefractionAngleCalc.h"
#include "CalcRefraction.h"
#include "TauConstant.h"
#include "SmoothRefrac.h"
#include "TruncRefracIn.h"
#include "TrimRefrac.h"
#include "RemoveMissingRefrac.h"
#include "GATS_Exception.h"
#include "TaperSolarShrink.h"
//#include "AvgExtentForRefrac.h"
#include "CalcRefracHigh.h"
#include "MergeRefrac.h" 
#include "Vec2Val.hpp"
#include "ValReversi.hpp"
#include "AddToEvent.hpp"

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



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



//
//-----------------------------------------------------------------------
// Utility Routines:
//-----------------------------------------------------------------------
//

using namespace std;


int RefractionAngleCalc(Event& L0, Event& L1, Event& Tmp, Event& SD, ConfigFile& cf)
{

    const string MODULE_NAME = "RefractionAngleCalc";

    const double missing = 2.0;
    RemoveMissingRefrac removeMissing(missing);
    
    valarray<double> refracTime;
    valarray<double> time;
    
    valarray<double> goodExtents;
    double exo_width;
    
//    AvgExtentForRefrac extentAveraging(100);
//    double meanExo;
    
    EventVar <double> timeArray;
    EventVar <double> sinkRateArray;
    EventVar <double> exoSolarExtent;
    EventVar <double> MeasuredExtent;
    EventVar <double> truncTime;
    EventVar <double> truncExtent;
    EventVar <double> eventVar;   
    
    valarray <double> refractionLow;
    valarray <double> refractionHigh;
    valarray <double> refractionMerged;

    TaperSolarShrink taper;

    double startTime;
    double stopTime;

    double rate; 
    double pixelConv;
    double t_0;
    double time_intv = 5.0;
    double timeRef;
    const double MISSING_VALUE = 10000;		// Large enough for the purpose of this module
    TruncRefracIn truncData;

    double beginEvent;
    double beginOcc;
    double endOcc;
    
    unsigned long indx;
    TauPointer tauFunctPtr;
    
    int retCode = 0;
    string strMsg;

    int taperSize = cf.GetInt(MODULE_NAME, "TaperLength");
    taper = TaperSolarShrink(taperSize);
    int smoothSize = cf.GetInt(MODULE_NAME, "SmoothingLength");

    double mergeCrit = cf.GetReal(MODULE_NAME, "MergeCriteria");
    int mergeHW = cf.GetInt(MODULE_NAME, "MergeHalfWidth");
        
    TrimRefrac trim(smoothSize);
    SmoothRefrac smooth(smoothSize);    

    MergeRefrac mergeFunc(mergeCrit, mergeHW);
    valarray <double> taperedShrinkage; 

    AddToEvent <double> addToL1Obj(L1);
    
    try {


        L1.getEventVar("L1_DetectorTimes", timeArray);	
        L1.getEventVar("L1_SolarExtent", exoSolarExtent); 
	
	    L1.getEventVar("SolarExtent", MeasuredExtent);
        L1.getEventVar("L1_SinkRate", sinkRateArray);
        rate = sinkRateArray[sinkRateArray.size()/2];


//      Only use data within occultation portion of measurement

        beginEvent = L1.getEventStartTime();
        beginOcc   = L1.getAtmosStartTime() - beginEvent; 
        endOcc     = L1.getAtmosStopTime() - beginEvent;
        TruncRefracIn truncData(beginOcc, endOcc);
        truncData(timeArray, MeasuredExtent, truncTime, truncExtent);
        
//	    Remove missing values
    	
	    removeMissing(truncTime, truncExtent, refracTime, goodExtents); 

         	
	    exo_width = exoSolarExtent[0];           
	    valarray<double> solar_ext(goodExtents.size()); 
        valarray<double>  smoothedExtents;
        valarray<double>  trimmedTime;
        valarray<double>  trimmedExtents;
        valarray<double>  solarShrinkage;        

//		solar_ext = exo_width*goodExtents;
        smooth(goodExtents, solar_ext);
        
         

//      taper input exo_width
         
        CalcRefraction calcRefr;
        double tau = fabs(exo_width/rate);
        tauFunctPtr.reset(new TauConstant(tau));  // Assume a constant sink rate... for now   
         
        calcRefr = CalcRefraction(tauFunctPtr);
        CalcRefracHigh highRefrac(tauFunctPtr);

        char modeFlag = L1.getSRSSFlag();

        cout << rate << endl;
        cout << exo_width << endl;
        
        if (modeFlag == 'r')
            {
            // Flip extent curve so that it looks like a sunset event.  Afterwards flip back.

            solar_ext = ValReversi(solar_ext);  
            trim(refracTime, solar_ext, trimmedTime, trimmedExtents);
            solarShrinkage.resize(trimmedExtents.size());
            refractionLow.resize(trimmedTime.size()); 
            refractionHigh.resize(trimmedTime.size());       
            refractionMerged.resize(trimmedTime.size());
            
//            meanExo = extentAveraging(trimmedExtents);
//            solarShrinkage = exo_width - trimmedExtents;
//            solarShrinkage = trimmedExtents.max() - trimmedExtents;
            solarShrinkage = exo_width - trimmedExtents;
            
            taper(solarShrinkage, taperedShrinkage);            
            refractionHigh = highRefrac.get_refraction(trimmedTime, solarShrinkage);            
            refractionLow = calcRefr.get_refraction(trimmedTime, taperedShrinkage);
            refractionMerged = mergeFunc.Combine(refractionHigh, refractionLow);
            
            refractionLow  = ValReversi(refractionLow);
            refractionHigh    = ValReversi(refractionHigh);
            taperedShrinkage = ValReversi(taperedShrinkage);  
            refractionMerged = ValReversi(refractionMerged);            
            }
         else
            {      
            trim(refracTime, solar_ext, trimmedTime, trimmedExtents);  
            solarShrinkage.resize(trimmedExtents.size());
            refractionLow.resize(trimmedTime.size());  
            refractionHigh.resize(trimmedTime.size());       
            refractionMerged.resize(trimmedTime.size());
                              
//            meanExo = extentAveraging(trimmedExtents);
//            solarShrinkage = exo_width - trimmedExtents;
//            solarShrinkage = trimmedExtents.max() - trimmedExtents;
            solarShrinkage = exo_width - trimmedExtents;

            taper(solarShrinkage, taperedShrinkage);            
            refractionHigh = highRefrac.get_refraction(trimmedTime, solarShrinkage);            
            refractionLow = calcRefr.get_refraction(trimmedTime, taperedShrinkage);
            refractionMerged = mergeFunc.Combine(refractionHigh, refractionLow);

            }

        eventVar = EventVar <double> ("", trimmedTime);
        eventVar.setName("L1_RefractionTime");
        addToL1Obj(eventVar);

        eventVar = EventVar <double> ("", taperedShrinkage);
        eventVar.setName("SolarShrinkage");
        addToL1Obj(eventVar);

        eventVar = EventVar <double> ("", refractionMerged);
        eventVar.setName("L1_RefractionAngles");
        addToL1Obj(eventVar);

        eventVar = EventVar <double> ("", refractionLow);
        eventVar.setName("RefractionAngles_LowAlgo");
        addToL1Obj(eventVar);

        eventVar = EventVar <double> ("", refractionHigh);
        eventVar.setName("RefractionAngles_HighAlgo");
        addToL1Obj(eventVar);
/*
std::cout << " in refractionanglecalc " << std::endl;
for(int i = 0; i< trimmedTime.size(); ++i) {
std::cout << refractionMerged[i] << "  " << trimmedTime[i] << std::endl; 
}
std::cout << "&" << std::endl;
*/
        return 0;
    
    } catch (exception &ex) {

        strMsg = string("Exception caught in ") + MODULE_NAME 
                        + string("SignalCorrection: ") + ex.what();
        L1.addLogEntry(strMsg);
        cerr << strMsg << endl;        
        retCode = -200;
                
    }
        
    return retCode;
}