//******************************************************************************
//** SCATMECH: Polarized Light Scattering C++ Class Library
//** 
//** File: grating.h
//**
//** Thomas A. Germer
//** Optical Technology Division, National Institute of Standards and Technology
//** 100 Bureau Dr. Stop 8443; Gaithersburg, MD 20899-8443
//** Phone: (301) 975-2876; FAX: (301) 975-6991
//** Email: thomas.germer@nist.gov
//**
//** Version: 6.00 (February 2008)
//**
//******************************************************************************
#ifndef SCATMECH_GRATING_H
#define SCATMECH_GRATING_H

#include "inherit.h"
#include "dielfunc.h"
#include <map>

namespace SCATMECH {
    
    class RCW_Model;
    
    class Grating : public Model {
    public:
        // Return the fourier component for specific order of a specific level.
        // level = 0 is the closest level to the incident direction.
        // recip == 0 returns the fourier component for epsilon.
        // recip == 1 returns the fourier component for 1/epsilon.
        virtual COMPLEX fourier(int order,int level,int recip=0) = 0;
                
        // Return the number of levels...
        virtual int get_levels()=0;
        
        // Return the thickness of a specific level (counting from TOP)...
        virtual double get_thickness(int level) = 0;
        
        void set_lambda(double _lambda) {lambda = _lambda; set_recalc(1);}
        double get_lambda() const {return lambda;}

        friend class RCW_Model;

    protected:
        
        // epsilon returns the dielectric function at lambda in form N-iK
        COMPLEX epsilon(const dielectric_function& e);
        
        double lambda;
        
    protected:

        DECLARE_MODEL();
        DECLARE_PARAMETER(double,period);
        DECLARE_PARAMETER(double,boundary);
        DECLARE_PARAMETER(dielectric_function,medium_i);
        DECLARE_PARAMETER(dielectric_function,medium_t);
    };
    
    void Register(const Grating* x);
    typedef Model_Ptr<Grating> Grating_Ptr;
    
    class Fourier_Component_Calculator {
    public:
        Fourier_Component_Calculator(double _period,int _order,int _recip,double _boundary);
        void enter(double x,COMPLEX y);
        COMPLEX result() { return sum; }
        
    private:
        double period,last_x,boundary;
        int order,recip;
        COMPLEX sum,last_y;
        bool begin;
    };
    
    class Generic_Grating : public Grating {
    public: 
        COMPLEX fourier(int order,int level,int recip=0);
        int get_levels() {SETUP(); return thickness.size();}
        double get_thickness(int level) {SETUP(); return thickness[level];}
                
        DECLARE_MODEL();
        DECLARE_PARAMETER(std::string,filename);
        DECLARE_PARAMETER(std::string,pstring);
        DECLARE_PARAMETER(int,nlayers);

    protected:
        void setup();
        
    private:
        std::vector<std::vector<std::string> > material;
        std::vector<std::vector<double> > position;
        std::vector<double> thickness;

        std::string filecontents;
        std::string old_filename;

        typedef std::map<std::string,COMPLEX> epsmap;
        epsmap eps;

        typedef std::map<std::string,double> varsmap;
        varsmap vars;
        
        std::vector<std::string> vnames;

        void error(const std::string& message) const;
        
        COMPLEX get_complex_value(std::istream& is) const;
    };
    
    class Null_Grating : public Grating {
        COMPLEX fourier(int i,int level,int recip=0) {return 0.;}
        int get_levels() {return 0;}
        double get_thickness(int level) {return 0;}
        DECLARE_MODEL();
    };

    class Single_Line_Grating : public Grating {
    public: 
        COMPLEX fourier(int order,int level,int recip=0);
        
    protected:
        void setup();
        int get_levels() {return nlevels;}
        double get_thickness(int level) {return height/nlevels;}
        
    private:
        COMPLEX eline;
        COMPLEX espace;
        
        DECLARE_MODEL();
        DECLARE_PARAMETER(dielectric_function,material); 
        DECLARE_PARAMETER(dielectric_function,space);
        DECLARE_PARAMETER(double,topwidth);
        DECLARE_PARAMETER(double,bottomwidth);
        DECLARE_PARAMETER(double,offset);
        DECLARE_PARAMETER(int,nlevels);
        DECLARE_PARAMETER(double,height);
    };
    
    class Corner_Rounded_Grating : public Grating {
    public:
        COMPLEX fourier(int order,int level,int recip=0);
        
        int get_levels() {
            return nlevels;
        }
        
        double get_thickness(int level) {
            SETUP();
            
            if (level==0) return height-(y[0]+y[1])/2.;
            if (level==nlevels-1) return (y[nlevels-2]+y[nlevels-1])/2.;
            if (level>0&&level<nlevels-1) return (y[level-1]+y[level])/2-(y[level]+y[level+1])/2.;
            return 0.;
        }
        
    protected:
        void setup();
        
    private:
        COMPLEX eps_line;
        
        DECLARE_MODEL();
        DECLARE_PARAMETER(double,width);
        DECLARE_PARAMETER(double,height);
        DECLARE_PARAMETER(dielectric_function,material);
        DECLARE_PARAMETER(double,sidewall);
        DECLARE_PARAMETER(double,radiusb);
        DECLARE_PARAMETER(double,radiust);
        DECLARE_PARAMETER(int,nlevels);
        
        std::vector<double> y;
        std::vector<double> x;
    };
    
    class Triangular_Grating : public Grating {
    public:
        COMPLEX fourier(int i,int level,int recip=0);
        
        int get_levels() {
            return nlevels;
        }
        
        double get_thickness(int level) {
            return amplitude/nlevels;
        }
        
    private:
        DECLARE_MODEL();
        DECLARE_PARAMETER(dielectric_function,material);
        DECLARE_PARAMETER(double,amplitude);
        DECLARE_PARAMETER(double,aspect);
        DECLARE_PARAMETER(int,nlevels);
    };
    
    class Sinusoidal_Grating : public Grating {
    public:
        COMPLEX fourier(int i,int level,int recip=0);
        
        int get_levels() {
            return nlevels;
        }
        
        double get_thickness(int level) {
            SETUP();
            int _nlevels = base!=0 ? nlevels-1 : nlevels;
            
            if (level==_nlevels) return base;
            
            double x1 = period/2.*level/_nlevels;
            double x0 = period/2.*(level+1)/_nlevels;
            return amplitude/2.*(cos(2*pi/period*x1)-cos(2*pi/period*x0));
        }
        
        DECLARE_MODEL();
        DECLARE_PARAMETER(double,amplitude);
        DECLARE_PARAMETER(double,base);
        DECLARE_PARAMETER(dielectric_function,material);
        DECLARE_PARAMETER(int,nlevels);
        
    };
    
}
#endif