DataClasses.h

Go to the documentation of this file.

/*
 *  SNOWPACK stand-alone
 *
 *  Copyright WSL Institute for Snow and Avalanche Research SLF, DAVOS, SWITZERLAND
*/
/*  This file is part of Snowpack.
    Snowpack is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Snowpack is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with Snowpack.  If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef DATACLASSES_H
#define DATACLASSES_H

#include <snowpack/SnowpackConfig.h>

#include <snowpack/Constants.h>
#include <meteoio/MeteoIO.h>

#include <string>
#include <vector>

enum {
        SOLID,  
        LIQUID, 
        GAS,   
        NUMBER_OF_PHASES
};

struct SNOW_OPTIC {
        double ggg;     
        double exteff;  
        double ssa;     
};

struct WL_STRUCT {
        double nm;      
        double perc;    
};

class ZwischenData {
        public:
                ZwischenData(): hoar24(48, 0.0), drift24(48, 0.0), hn3(144, 0.0), hn24(144, 0.0) {}
                void reset();                

                friend std::ostream& operator<<(std::ostream& os, const ZwischenData& data);
                friend std::istream& operator>>(std::istream& is, ZwischenData& data);

                std::vector<double> hoar24;  
                std::vector<double> drift24; 
                std::vector<double> hn3;    
                std::vector<double> hn24;   
};

class CurrentMeteo {
        public:
                CurrentMeteo();
                CurrentMeteo(const SnowpackConfig& i_cfg);
                
                void reset(const SnowpackConfig& i_cfg);
                void setMeasTempParameters(const mio::MeteoData& md);
                size_t getNumberMeasTemperatures() const;
                size_t getNumberFixedRates() const;
                size_t getMaxNumberMeasTemperatures() const;
                void getFixedPositions(std::vector<double>& positions) const;
                size_t getNumberFixedPositions() const;
                void copySnowTemperatures(const mio::MeteoData& md, const unsigned int& current_slope);
                void copySolutes(const mio::MeteoData& md, const size_t& i_number_of_solutes);

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const CurrentMeteo& data);
                friend std::istream& operator>>(std::istream& is, CurrentMeteo& data);

                mio::Date date;  
                double ta;       
                double rh;       
                double rh_avg;   
                double vw;       
                double vw_avg;   
                double vw_max;   
                double dw;       
                double vw_drift; 
                double dw_drift; 
                double ustar;    
                double z0;       
                double psi_s;    
                double iswr;     
                double rswr;     
                double mAlbedo;  
                double diff;     
                double dir_h;    
                double elev;     
                double ea;       
                double tss;      
                double tss_a12h; 
                double tss_a24h; 
                double ts0;      
                double psum;     
                double psum_ph;  
                double psum_tech;
                double hs;       
                double hs_a3h;   
                double hs_rate;  
                double adv_heat; 

                std::vector<double> ts;    
                std::vector<double> zv_ts; 
                std::vector<double> conc;  
                double rho_hn;             

        private:
                size_t getNumberMeasTemperatures(const mio::MeteoData& md);

                std::vector<double> fixedPositions; 
                double minDepthSubsurf;             
                size_t maxNumberMeasTemperatures;   
                size_t numberMeasTemperatures;      
                size_t numberFixedRates;
};

enum SN_FIELD{
        TEMPERATURE, 
        SEEPAGE,     
        SETTLEMENT,  
        N_SN_FIELDS
};

enum {
        SOIL,  
        ICE,   
        WATER, 
        AIR,   
        N_COMPONENTS
};

enum SN_SOIL_DATA{
        SOIL_RHO, 
        SOIL_K,   
        SOIL_C,   
        N_SOIL_FIELDS
};

class LayerData {
        public:
                LayerData();

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const LayerData& data);
                friend std::istream& operator>>(std::istream& is, LayerData& data);

                mio::Date depositionDate;   
                double hl;                  
                size_t ne;                  
                double tl;                  
                double phiSoil;             
                double phiIce;              
                double phiWater;            
                double phiVoids;            
                std::vector<double> cSoil;  
                std::vector<double> cIce;   
                std::vector<double> cWater; 
                std::vector<double> cVoids; 
                double SoilRho;             
                double SoilK;               
                double SoilC;               
                double rg;                  
                double sp;                  
                double dd;                  
                double rb;                  
                unsigned short int mk;      
                double hr;                  
                double CDot;                
                double metamo;              
};

class SN_SNOWSOIL_DATA {
        public:
                SN_SNOWSOIL_DATA() : meta(), profileDate(), nN(0), Height(0.),
                     nLayers(0), Ldata(), HS_last(0.), Albedo(0.), SoilAlb(0.), BareSoil_z0(0.),
                     Canopy_Height(0.), Canopy_LAI(0.), Canopy_BasalArea(0.004), Canopy_Direct_Throughfall(0.),
                     WindScalingFactor(1.), ErosionLevel(0), TimeCountDeltaHS(0.) {}

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const SN_SNOWSOIL_DATA& data);
                friend std::istream& operator>>(std::istream& is, SN_SNOWSOIL_DATA& data);

                mio::StationData meta;            
                mio::Date profileDate;            
                size_t nN;                        
                double Height;                    
                size_t nLayers;                   
                std::vector<LayerData> Ldata;     
                double HS_last;                   
                double Albedo;                    
                double SoilAlb;                   
                double BareSoil_z0;               
                double Canopy_Height;             
                double Canopy_LAI;                
                double Canopy_BasalArea;          
                double Canopy_Direct_Throughfall; 
                double WindScalingFactor;         
                int    ErosionLevel;              
                double TimeCountDeltaHS;          
};

class ElementData {
        public:
                typedef enum YOUNG_MODULUS {
                            Sigrist, 
                            Pow, 
                            Exp 
                } Young_Modulus;
                
                ElementData(const unsigned short int& in_ID);
                
                bool checkVolContent() const;
                void heatCapacity();
                double coldContent() const;
                double extinction() const;
                void opticalEquivalentGrainSize();
                void snowResidualWaterContent();
                static double snowResidualWaterContent(const double& theta_i);
                double soilFieldCapacity() const;

                double snowElasticity() const;
                double neckStressEnhancement() const;
                double concaveNeckRadius() const;
                double neckLength() const;
                double neck2VolumetricStrain() const;

                void snowType();
                unsigned short int getSnowType() const;
                static unsigned short int snowType(const double& dendricity, const double& sphericity, const double& grain_dia, const unsigned short int& marker,
                        const double& theta_w, const double& res_wat_cont);
                static double getYoungModule(const double& rho_slab, const Young_Modulus& model);

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const ElementData& data);
                friend std::istream& operator>>(std::istream& is, ElementData& data);

                mio::Date depositionDate;  
                double L0, L;              
                double Te;                 
                double gradT;              
                double melting_tk;         
                double freezing_tk;        
                std::vector<double> theta; 
                mio::Array2D<double> conc; 
                std::vector<double> k;     
                //   Stored in order to visualize constitutive laws
                //   Will be used for creep field hydraulic conductivity in m3 s kg-1
                std::vector<double> c;     
                //   Will also be used for creep specific snow water capacity  in m3 J-1
                std::vector<double> soil;  
                double Rho;                
                double M;                  
                double sw_abs;             
                // Snow Metamorphism Data
                double rg;                 
                double dd;                 
                double sp;                 
                double ogs;                
                double rb;                 
                double N3;                 
                size_t mk;                 
                unsigned short int type;   
                double metamo;             
                double dth_w;              
                double res_wat_cont;       
                double Qmf;                
                double QIntmf;             
                double dEps, Eps, Eps_e, Eps_v;      
                double Eps_Dot, Eps_vDot;        
                double E;                  
                double S;                  
                double C;                  
                double CDot;               
                double ps2rb;              
                double s_strength;         
                double hard;               
                double S_dr;               
                double crit_cut_length;    
                double theta_r;            
                double lwc_source;         
                double SlopeParFlux;       
                //NIED (H. Hirashima)
                double dhf;
                
                unsigned short int ID;    
                static const unsigned short int noID;
};

class NodeData {
        public:
                NodeData() : z(0.), u(0.), f(0.), udot(0.), T(0.), S_n(0.), S_s(0.), ssi(6.), hoar(0.),
                             dhf(0.), S_dhf(0.), Sigdhf(0.) {} //HACK: set ssi to max_stability!

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const NodeData& data);
                friend std::istream& operator>>(std::istream& is, NodeData& data);

                double z;    
                double u;    
                double f;    
                double udot; 
                double T;    
                double S_n;  
                double S_s;  
                double ssi;  
                double hoar; 

                //NIED (H. Hirashima)
                double dhf;
                double S_dhf;
                double Sigdhf;
};

class CanopyData {
        public:
                CanopyData() : storage(0.), temp(0.), sigf(0.), ec(0.), lai(0.), z0m(0.), z0h(0.), zdispl(0.),
                     height(0.), direct_throughfall(0.), ra(0.), rc(0.), rs(0.), rstransp(0.), canopyalb(0.),
                     totalalb(0.), wetfraction(0.), intcapacity(0.), rswrac(0.), iswrac(0.), rswrbc(0.),
                     iswrbc(0.), ilwrac(0.), rlwrac(0.), ilwrbc(0.), rlwrbc(0.), rsnet(0.), rlnet(0.),
                     sensible(0.), latent(0.), latentcorr(0.), transp(0.), intevap(0.),
                     interception(0.), throughfall(0.), snowunload(0.),
                     snowfac(0.), rainfac(0.),liquidfraction(0.),
                     sigftrunk(0), Ttrunk(0.), CondFluxCanop(0.), CondFluxTrunks(0.),
                     LWnet_Trunks(0.), SWnet_Trunks(0.), QStrunks(0.),
                     forestfloor_alb(0.), BasalArea(0), HMLeaves(0.), HMTrunks(0.) {}

                void reset(const bool& cumsum_mass);
                void initializeSurfaceExchangeData();
                void multiplyFluxes(const double& factor);

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const CanopyData& data);
                friend std::istream& operator>>(std::istream& is, CanopyData& data);

                // State variable
                double storage;     
                double temp;        
                double sigf;        
                double ec;          
                // parameters
                double lai;
                double z0m;
                double z0h;
                double zdispl;
                double height;
                double direct_throughfall;
                // aerodynamic resistances
                double ra;          
                double rc;          
                double rs;          
                double rstransp;    
                // Averaged variables
                double canopyalb;   
                double totalalb;    
                double wetfraction; 
                double intcapacity; 
                // Radiations
                double rswrac;      
                double iswrac;      
                double rswrbc;      
                double iswrbc;      
                double ilwrac;      
                double rlwrac;      
                double ilwrbc;      
                double rlwrbc;      
                double rsnet;       
                double rlnet;       
                // Turbulent fluxes
                double sensible;
                double latent;
                double latentcorr;
                // Evap fluxes
                double transp;
                double intevap;
                // Mass fluxes
                double interception;
                double throughfall;
                double snowunload;
                
                double snowfac;     
                double rainfac;     
                double liquidfraction;
                double sigftrunk;   
                double Ttrunk;      
                double CondFluxCanop; 
                double CondFluxTrunks; 
                double LWnet_Trunks; 
                double SWnet_Trunks; 
                double QStrunks;      
                double forestfloor_alb; 
                double BasalArea;    
                double HMLeaves;     
                double HMTrunks;     

};

class SnowStation {
        public:
                explicit SnowStation(const bool& i_useCanopyModel=true, const bool& i_useSoilLayers=true);
                SnowStation(const SnowStation& c);

                ~SnowStation();
                SnowStation& operator=(const SnowStation&); 

                void initialize(const SN_SNOWSOIL_DATA& SSdata, const size_t& i_sector);
                void resize(const size_t& number_of_elements);

                void reduceNumberOfElements(const size_t& rnE);
                void combineElements(const size_t& number_top_elements, const bool& reduce_n_elements);
                void combineElements(const size_t& number_top_elements, const bool& reduce_n_elements, const size_t& cond);
                static bool combineCondition(const ElementData& Edata0, const ElementData& Edata1, const double& depth, const bool& reduce_n_elements);
                static void mergeElements(ElementData& Edata0, const ElementData& Edata1, const bool& merge, const bool& topElement);
                void splitElement(const size_t& e);     //Split an element
                void splitElements();                   //Check for splitting, calls splitElement(...) for actual splitting

                void compSnowpackMasses();
                void compSnowpackInternalEnergyChange(const double& sn_dt);
                void compSoilInternalEnergyChange(const double& sn_dt);
                double getLiquidWaterIndex() const;
                double getModelledTemperature(const double& z) const;

                size_t getNumberOfElements() const;
                size_t getNumberOfNodes() const;
                bool isGlacier(const bool& hydro=false) const;
                bool hasSoilLayers() const;

                size_t find_tag(const size_t& tag) const;

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const SnowStation& data);
                friend std::istream& operator>>(std::istream& is, SnowStation& data);

                mio::StationData meta;      
                double cos_sl;              
                size_t sector;              

                CanopyData Cdata;           
                static const void* Seaice;  
                double pAlbedo;             
                double Albedo;              
                double SoilAlb;             
                double BareSoil_z0;         
                size_t SoilNode;            
                double Ground;              
                double cH;                  
                double mH;                  
                double mass_sum;            
                double swe;                 
                double lwc_sum;             
                double hn;                  
                double rho_hn;              
                size_t ErosionLevel;        
                double ErosionMass;         
                char S_class1;               
                char S_class2;               
                double S_d;                 
                double z_S_d;               
                double S_n;                 
                double z_S_n;               
                double S_s;                 
                double z_S_s;               
                double S_4;                 
                double z_S_4;               
                double S_5;                 
                double z_S_5;               
                std::vector<NodeData> Ndata;    
                std::vector<ElementData> Edata; 
                void *Kt;                   
                double ColdContent;         
                double ColdContentSoil;     
                double dIntEnergy;          
                double dIntEnergySoil;      
                double meltFreezeEnergy;    
                double meltFreezeEnergySoil;
                double ReSolver_dt;         
                bool windward;              
                double WindScalingFactor;   
                double TimeCountDeltaHS;    

                static const double comb_thresh_l, comb_thresh_ice, comb_thresh_water;
                static const double comb_thresh_dd, comb_thresh_sp, comb_thresh_rg;
                static const double thresh_moist_snow, thresh_moist_soil;
                static const size_t number_top_elements;
                static unsigned short number_of_solutes;  

        private:
                static double flexibleMaxElemLength(const double& depth); 
                
                size_t nNodes;                      
                size_t nElems;                      
                unsigned short int maxElementID;    
                bool useCanopyModel, useSoilLayers; 
};

class BoundCond {

        public:
                BoundCond() : lw_out(0.), lw_net(0.), qs(0.), ql(0.), qr(0.), qg(Constants::undefined) {}
                const std::string toString() const;
                void reset();

                double lw_out;  
                double lw_net;  
                double qs;      
                double ql;      
                double qr;      
                double qg;      
};

class SurfaceFluxes {
        public:
                enum SN_MASS_CHANGES {
                        MS_TOTALMASS,      
                        MS_SWE,            
                        MS_WATER,          
                        MS_HNW,            
                        MS_RAIN,           
                        MS_WIND,           
                        MS_EVAPORATION,    
                        MS_SUBLIMATION,    
                        MS_SNOWPACK_RUNOFF,
                        MS_SOIL_RUNOFF,    
                        N_MASS_CHANGES     
                };

                const std::string toString() const;
                friend std::ostream& operator<<(std::ostream& os, const SurfaceFluxes& data);
                friend std::istream& operator>>(std::istream& is, SurfaceFluxes& data);

                SurfaceFluxes();

                void reset(const bool& cumsum_mass);
                void compSnowSoilHeatFlux(const SnowStation& Xdata);
                void collectSurfaceFluxes(const BoundCond& Bdata, SnowStation& Xdata, const CurrentMeteo& Mdata);
                void multiplyFluxes(const double& factor);

                double lw_in;      
                double lw_out;     
                double lw_net;     
                double qs;         
                double ql;         
                double hoar;       
                double qr;         
                double qg;         
                double qg0;        
                double sw_hor;     
                double sw_in;      
                double sw_out;     
                double qw;         
                double sw_dir;     
                double sw_diff;    
                double pAlbedo;    
                double mAlbedo;    
                double dIntEnergy;           
                double dIntEnergySoil;       
                double meltFreezeEnergy;     
                double meltFreezeEnergySoil; 

                double drift;      
                std::vector<double> mass; 
                std::vector<double> load; 
                double dhs_corr;   
                double cRho_hn;    
                double mRho_hn;    

};

//HACK: could it be moved to plugins? (as well as Aggregate)
class SnowProfileLayer {
        public:
                SnowProfileLayer();

                void average(const double& w1, const double& w2, const SnowProfileLayer& Pdata);
                static std::vector<SnowProfileLayer> generateProfile(const mio::Date& dateOfProfile, const SnowStation& Xdata, const double hoar_density_surf, const double hoar_min_size_surf);

                // Profile meta data
                mio::Date profileDate; 
                std::string stationname;
                unsigned char  loc_for_snow;
                unsigned char  loc_for_wind;

                mio::Date depositionDate;   
                double height;         
                double rho;            
                double T;              
                double gradT;          
                double v_strain_rate;  
                double theta_i;        
                double theta_w;        
                double theta_a;        
                double grain_size;     
                double bond_size;      
                double dendricity;     
                double sphericity;     
                double ogs;            
                double coordin_num;    
                unsigned short int marker;         
                short unsigned int type; 
                double hard;           

        private:
                void generateLayer(const ElementData& Edata, const NodeData& Ndata);
                void generateLayer(const ElementData& Edata, const NodeData& Ndata,
                                   const mio::Date& dateOfProfile, const double hoar_density_surf);
};

class RunInfo {
        public:
                RunInfo();
                RunInfo(const RunInfo& orig);
                RunInfo& operator=(const RunInfo&) {return *this;} //everything is static, so we can not change anything

                const std::string version;   
                const mio::Date computation_date; 
                const std::string compilation_date; 
                const std::string user; 
                const std::string hostname; 

        private:
                static mio::Date getRunDate();
                static std::string getCompilationDate();
};

#endif