using System;
using System.Collections.Generic;
using System.Text;
using System.IO;

namespace GeoFly
{
    /// <summary>
    /// 天气发生器类
    /// </summary>
    public class ClimateGene
    {
        Random Random = new Random();
        public ClimateGene()
        {
        }
        /// <summary>
        /// 生成气温数据
        /// </summary>
        /// <param name="datetype"></param>
        /// <param name="simulatecode"></param>
        /// <returns></returns>
        public void TempSimulate(int datetype, int simulatecode)
        {
            ////this.CalcPcpCoFactor(pcpcode);

            //List<string> saArray;
            //int strnum;

            ////从到1天到总天数，开始模拟
            //for (int i = 0; i < totdy; i++)
            //{
            //    //获取第i天的日期时间
            //    //DateTime strdate = this.PcpData.wgnDate[i];
            //    //获取月
            //    //int imon = strdate.Month;

            //    //this.TempResiduals();
            //    //DailyPcp(imon, strdate);
            //    //DailyTmp(imon, strdate);				//simulate daily max/min temperature
            //    //DailySlr(imon, strdate);				//simulate daily solar radiation
            //    //DailyWnd(imon, strdate);				//simulate daily wind 
            //    //DailyHmd(imon, strdate);				//simulate daily relative humidity
                
            //}           
        }     
        /// <summary>
        ///使用偏正态分布模型(见SWAT，刘永和写)：generates precipitation data when the user chooses to simulate 
        ///or when data is missing for particular days in the weather file
        /// </summary>
        /// <param name="strdate">日期</param>
        /// <param name="IsDryPreDay">前一天是否为干天的标志</param>
        /// <returns>当前日期的降水量</returns>
        public double GenPcp_Skewed(DateTime strdate, bool IsDryPreDay)
        {
            double rand1 = Random.NextDouble();
            int month = strdate.Month-1;
            if (IsDryPreDay) //前一天无雨，则执行PDW概率
            {
                if (rand1 <= this.pData[month].PDW)  //小于这个概率为雨天
                {
                    double delta = this.pData[month].PCPSTD; //标准差
                    double gmon = this.pData[month].PCPSKW;  //偏态系数
                    double miu = this.pData[month].MPDM;
                    double SNDday = this.GaussRandom();
                    double nomi = Math.Pow(((SNDday - gmon / 6.0) * gmon / 6.0 + 1), 3) - 1;
                    double result= miu + 2.0 * delta * nomi / gmon;
                    if (result < 0.1)
                        result = 0.1;
                    return result;
                }
                else
                {
                    return 0;  //否则无雨
                }
            }
            else //前一天为雨天，则执行PWW概率
            {
                if (rand1 <= this.pData[month].PWW)  //小于这个概率为雨天
                {
                    double delta = this.pData[month].PCPSTD; //标准差
                    double gmon = this.pData[month].PCPSKW;  //偏态系数
                    double miu = this.pData[month].MPDM;
                    double SNDday = this.GaussRandom();
                    double nomi = Math.Pow(((SNDday - gmon / 6.0) * gmon / 6.0 + 1), 3) - 1;
                    double result= miu + 2.0 * delta * nomi / gmon;
                    if (result < 0.1)
                        result = 0.1;
                    return result;
                }
                else
                {
                    return 0;  //否则无雨
                }
            }
        }
        /// <summary>
        /// 采用Richardson提到的指数分布来生成日降水量，lamda用遗传算法求出
        /// </summary>
        /// <param name="strdate"></param>
        /// <param name="IsDryPreDay"></param>
        /// <returns></returns>
        public double GenPcp_expon(DateTime strdate, bool IsDryPreDay,double lamda)
        {
            throw new Exception("");
        }
        /// <summary>
        /// 根据SWAT的天气发生器计算气温(刘永和修改)
        /// 主要包括日最高气温、日平均气温
        /// </summary>
        /// <param name="mon"></param>
        /// <param name="strdate"></param>
        /// <param name="WetDays">月内的雨天天数(根据生成的雨量数据可以先统计出来)</param>
        /// <param name="TotalDays">月的总天数</param>
        /// <param name="bt">scaling factor(SWAT)</param>
        /// <returns>返回日内均温和日内最高温。均温下标为0，最高温的下标为1</returns>
        public double[] TempGen(int mon, DateTime strdate,bool IsDryDay, int WetDays, int TotalDays,double bt)
        {
            Random rnd = new Random();       
            
            double TMPMX = this.pData[mon].TMPMX;            
            double TMPMN = this.pData[mon].TMPMN;            
            double TMPSTDMN = this.pData[mon].TMPSTDMN;
            double TMPSTDMX = this.pData[mon].TMPSTDMX;
            double PR_W = this.pData[mon].PR_W;
            //计算干天日平均最高温度
            double uDmx = TMPMX + bt * ((double)WetDays / TotalDays) * (TMPMX - TMPMN);
            //计算雨天日平均最高温度
            double uWmx=uDmx-bt*(TMPMX - TMPMN);
            double[] dayTemp = new double[2];
            dayTemp[0] = TMPMN + wgncur[1] * TMPSTDMN; //生成的日内均温
            
            if (IsDryDay)
            {
                dayTemp[1] = uDmx + wgncur[0] * TMPSTDMX;
            }
            else
            {
                dayTemp[1] = uWmx + wgncur[0] * TMPSTDMX; //生成的日内最高温
            }
            return dayTemp;
        }
        /// <summary>
        /// 生成Solor radiation
        /// </summary>
        /// <param name="mon"></param>
        /// <param name="strdate"></param>
        /// <param name="dscale">默认为0.8</param>
        public double SlrGen(int mon, DateTime strdate,bool IsDryDay,int WetDays, int TotalDays, double br,double dscale,double lat)
        {
            //月内平均日太阳辐射
            double urad=this.pData[mon].SOLARAV;
            //月内干天平均日太阳辐射
            double uDrad = urad * TotalDays / (br * WetDays + (TotalDays - WetDays));
            //月内雨天平均日太阳辐射
            double uWrad = br * uDrad;
            //下面计算日最大太阳辐射
            int dn = DateAndTime.GetDnInYear(strdate.Year, strdate.Month, strdate.Day);
            SolarRadiation slrrad = new SolarRadiation(strdate.Year, dn, lat);
            double MaxDailySolar=slrrad.ExtraTerrRad(dn,lat)*dscale;
            //计算月内日太阳辐射的标准差
            double delta_rad=(MaxDailySolar-urad)/4;
            double dayslr = -9999;
            if (IsDryDay)
            {
                dayslr = uDrad + wgncur[2] * delta_rad;
            }
            else
            {
                dayslr = uWrad + wgncur[2] * delta_rad;
            }
            return dayslr;
        }
        /// <summary>
        /// 生成strdate那天的风速数据(刘永和写)
        /// </summary>
        /// <param name="mon"></param>
        /// <param name="strdate"></param>
        public double WndGen(int mon, DateTime strdate)
        {
            Random rnd = new Random();            
            return this.pData[mon].WNDAV * Math.Pow(-Math.Log(rnd.NextDouble()), 0.3);
        }
        /// <summary>
        /// 生成相对湿度(刘永和写)
        /// </summary>
        /// <param name="mon"></param>
        /// <param name="strdate"></param>
        public double HmdGen(int mon, DateTime strdate)
        {           
            
            //月平均气温
            double TmpMeanMon = (this.pData[mon].TMPMX + this.pData[mon].TMPMN) / 2.0;
            //月平均露点温度
            double DEWPT = this.pData[mon].DEWPT;
            //计算露点温度下的饱和水汽压，即实际水汽压(当前未必是饱和的)
            double Emon = Math.Exp((16.78 * DEWPT - 116.9) / (DEWPT + 237.3));
            //计算月平均饱和水气压
            double Emean = Math.Exp((16.78 * TmpMeanMon - 116.9) / (TmpMeanMon + 237.3));
            //计算月平均相对湿度：实际水汽压/饱和水汽压
            double Rhmon = Emon / Emean;
            //相对湿度最大上限
            double RhmonU = Rhmon + (1 - Rhmon) * Math.Exp(Rhmon - 1);
            //相对湿度最小下限
            double RhmonL = Rhmon * (1 - Math.Exp(-Rhmon));
            //概率阈值
            double threshold = (Rhmon - RhmonL) / (RhmonU - RhmonL);
            //下面套用SWAT的三角分布公式来计算日相对湿度
            double Rhmon_mean=(RhmonL+Rhmon+RhmonU)/3;
            double rnd1 = Random.NextDouble();
            double dayhmd = -9999;
            if (rnd1 <= threshold)
            {
                 dayhmd = Rhmon * (RhmonL + Math.Sqrt(rnd1 * (RhmonU - RhmonL) * (Rhmon - RhmonL))) / Rhmon_mean;
            }
            else
            {
                double sqrtValue=Math.Sqrt((RhmonU-RhmonL)*(1-rnd1)/(RhmonU-Rhmon));
                 dayhmd=Rhmon*(RhmonU-(RhmonU-Rhmon)*sqrtValue/Rhmon_mean);
            }
            return dayhmd;
        }
        /// <summary>
        /// 计算残差：RS(i)=A*RS(i-1)+B*E(i)
        /// </summary>
        public void TempResiduals()
        {
            double[] xx = new double[3];
            for (int i = 0; i < 3; i++)
            {
                wgnold[i] = wgncur[i]; //迭代式赋值，初始化
                xx[i] = 0;
                e[i] = 0;
            }            

            e[0] = GaussRandom();
            e[1] = GaussRandom();
            e[2] = GaussRandom();

            for (int i = 0; i < 3; i++)
            {
                wgncur[i] = 0;
                for (int j = 0; j < 3; j++)
                {
                    wgncur[i] += e[j] * MatrixFuncs.B[i, j];
                    xx[i] += wgnold[j] * MatrixFuncs.A[i, j];
                }
            }
            for (int i = 0; i < 3; i++)
            {
                wgncur[i] += xx[i];
            }
        }
        /// <summary>
        /// 计算日最高气温的月均值的多年再平均(度)
        /// </summary>
        private void CalTempMax()
        {
            if (this.LoadMaxTemprature == false)
                throw new Exception("还没有导入最高气温数据");
            int YearsCount =this.WeatherData.Keys.Count;
            //下面统计日最高气温的月均值的多年平均
            //--------------------------------------------------------
            double[] YearlyTotal = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlyTotal[month] += this.GetMean(WeatherData[year].MaxTemperature, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].TMPMX = YearlyTotal[month] / YearsCount;
            }
        }
        /// <summary>
        /// 日最高气温的月均值的多年再平均(度)
        /// </summary>
        private void CalTempMean()
        {
            if (this.LoadMeanTemprature == false)
                throw new Exception("还没有导入平均温度数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日平均气温的月均值的多年平均
            //----------------------------------------------------------------------
            double[] YearlyTotal = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlyTotal[month] += this.GetMean(WeatherData[year].MeanTemperature, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].TMPMN = YearlyTotal[month] / YearsCount;
            }
        }
        /// <summary>
        /// 日太阳辐射的月均值的多年平均
        /// </summary>
        private void CalSolarRadMean()
        {
            if (this.LoadSolarRadiation == false)
                throw new Exception("还没有导入太阳辐射数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日平均气温的月均值的多年平均
            double[] YearlySum = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlySum[month]+= this.GetMean(WeatherData[year].SolarRadiation, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].SOLARAV = YearlySum[month] / YearsCount;
            }
        }
        /// <summary>
        /// 日风速的月均值的多年平均
        /// </summary>
        private void CalWindSpeedMean()
        {
            if (this.LoadWindSpeed == false)
                throw new Exception("还没有导入风速数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            //----------------------------------------------------------------------
            double[] YearlyTotal = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlyTotal[month] += this.GetMean(WeatherData[year].WindSpeed, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].WNDAV = YearlyTotal[month] / YearsCount;
            }
            
        }
        /// <summary>
        /// 日降水量的月均值的多年平均
        /// </summary>
        private void CalPreciMean()
        {
            if (this.LoadPrecipitation == false)
                throw new Exception("还没有导入降水数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            for (int month = 0; month < 12; month++)
            {
                double YearlyTotal =0;
                foreach (int year in WeatherData.Keys)
                {
                    YearlyTotal += this.GetMean(WeatherData[year].Precipitation, month);
                }
                //日降水量的月均值的多年再平均
                this.pData[month].MPDM = YearlyTotal / YearsCount;
            }
        }
        /// <summary>
        /// 获取month月数据的标准差
        /// </summary>
        /// <param name="data"></param>
        /// <param name="month"></param>
        /// <returns></returns>
        public double GetSTD(double[,] data,int month)
        {
            double mean = this.GetMean(data, month);
            double std = 0;
            for (int day = 0; day < daysInMonth[month]; day++)
            {
                double value = data[month, day] - mean;
                std += value * value;
            }
            std /= daysInMonth[month];
            return Math.Sqrt(std);
        }
        /// <summary>
        /// 获取month月的平均值
        /// </summary>
        /// <param name="data"></param>
        /// <param name="month"></param>
        /// <returns></returns>
        public double GetMean(double[,] data, int month)
        {
            double sum = 0;
            for (int day = 0; day < daysInMonth[month]; day++)
            {
                sum += data[month, day];
            }
            return sum /daysInMonth[month];
        }
        public int[] daysInMonth = { 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31 };
        /// <summary>
        /// 日最高温度的月方差的多年平均值
        /// </summary>
        private void CalTempMaxSTD()
        {
            if (this.LoadMaxTemprature == false)
                throw new Exception("还没有导入最高温度数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            //----------------------------------------------------------------------
            double[] YearlySum = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlySum[month] += this.GetSTD(WeatherData[year].MaxTemperature, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].TMPSTDMX = YearlySum[month] / YearsCount;
            }
        }
        /// <summary>
        /// 日平均温度的月方差的多年平均值
        /// </summary>
        private void CalTempMeanSTD()
        {
            if (this.LoadMeanTemprature == false)
                throw new Exception("还没有导入平均温度数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            //----------------------------------------------------------------------
            double[] YearlySum = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlySum[month] += this.GetSTD(WeatherData[year].MeanTemperature, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].TMPSTDMN = YearlySum[month] / YearsCount;
            }
        }
        /// <summary>
        /// 日降水量的月方差的多年平均值
        /// </summary>
        private void CalPreciSTD()
        {
            if (this.LoadPrecipitation == false)
                throw new Exception("还没有导入降水数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            
            for (int month = 0; month < 12; month++)
            {
                double YearlySum = 0;
                foreach (int year in WeatherData.Keys)
                {
                    YearlySum += this.GetSTD(WeatherData[year].Precipitation, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].PCPSTD = YearlySum/ YearsCount;
            }
        }
        /// <summary>
        /// 日露点温度的月方差的多年平均值
        /// </summary>
        private void CalDewPoint()
        {
            if (this.LoadDewPoint == false)
                throw new Exception("还没有计算露点温度数据");
            int YearsCount = WeatherData.Keys.Count;
            //下面统计日风速的月均值的多年平均
            double[] YearlySum = new double[12];
            for (int month = 0; month < 12; month++)
            {
                foreach (int year in WeatherData.Keys)
                {
                    YearlySum[month] += this.GetSTD(WeatherData[year].DewPoint, month);
                }
                //日最高气温的月均值的多年再平均
                this.pData[month].DEWPT = YearlySum[month] / YearsCount;
            }
        }
        /// <summary>
        /// 计算降水偏态系数
        /// </summary>
        private void CalPCPSkew()
        {
            if (this.LoadPrecipitation == false)
                throw new Exception("还没有导入降水数据");
            int YearsCount = WeatherData.Keys.Count;
            
            for (int month = 0; month < 12; month++)
            {
                double YearlySum = 0;
                foreach (int year in WeatherData.Keys)
                {
                    double mean= this.GetMean(WeatherData[year].Precipitation, month);
                    double std = this.GetSTD(WeatherData[year].Precipitation, month);
                    double Cv = std / mean;//计算变差系数，即标准差与均值的比值
                    double sum = 0;
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        double dist = WeatherData[year].Precipitation[month, day] - mean;
                        sum+=Math.Pow(dist, 3);
                    }
                    
                    //sum *= daysInMonth[month];
                    //double denomi=(daysInMonth[month] - 1) * (daysInMonth[month] - 2) * Math.Pow(std, 3);
                    if(sum!=0) //当sum为0时denomi必然为0；当sum不为0时，denomi也不为0
                        YearlySum += sum / this.daysInMonth[month] / Math.Pow(mean, 3) / Math.Pow(Cv, 3);
                }
                if (YearsCount == 0)
                    throw new Exception("");
                this.pData[month].PCPSKW = YearlySum/ YearsCount;
            }
        }
        /// <summary>
        /// 统计多年月平均降水天数
        /// </summary>
        private void CalPCPDays()
        {
            if (this.LoadPrecipitation == false)
                throw new Exception("还没有导入降水数据");
            int YearsCount = WeatherData.Keys.Count;            
            for (int month = 0; month < 12; month++)
            {
                int total_wetdays=0;
                foreach (int year in WeatherData.Keys)
                {
                    int count=0;
                    for(int day=0;day<this.daysInMonth[month];day++)
                    {
                        if(WeatherData[year].Precipitation[month,day]>0.1)
                            count+=1;
                    }
                    total_wetdays+=count;
                }
                this.pData[month].PCPD = total_wetdays / YearsCount;
            }
        }
        /// <summary>
        /// 统计多年月总降水量的均值
        /// </summary>
        private void CalMeanPCPMonSum()
        {
            int YearsCount = WeatherData.Keys.Count;
            for (int month = 0; month < 12; month++)
            {
                double All = 0;
                foreach (int year in WeatherData.Keys)
                {
                    double MonthlySum = 0;//存放月总降水量
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        MonthlySum+=WeatherData[year].Precipitation[month,day];
                    }
                    All += MonthlySum;
                }
                this.pData[month].PCPMM = All / YearsCount;
            }
        }
        /// <summary>
        /// 统计前一天为干天当前天为雨天的概率
        /// </summary>
        private void CalProDW()
        {
            
            for (int month = 0; month < 12; month++)
            {
                int PreMonth = month - 1;
                double All = 0;
                int YearsCount = 0;
                foreach (int year in WeatherData.Keys)
                {
                    int DryDays = 0;//这个月的干天数
                    bool IsPreDry = false;
                    if (PreMonth == -1) //如果数据中前一个月没法对证，就认为IsPreDry=true;
                    {
                        IsPreDry = true;
                    }
                    else if (WeatherData[year].Precipitation[PreMonth, this.daysInMonth[PreMonth] - 1] < 0.1)  //从上个月的最后一天开始
                    {
                        IsPreDry = true;
                    }
                    int DWCount = 0;
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        if (IsPreDry) //统计前一天为干天后一天为雨天的次数
                        {
                            DryDays += 1;
                            if (WeatherData[year].Precipitation[month, day] >= 0.1)
                                DWCount += 1;
                        }
                        if (WeatherData[year].Precipitation[month, day] >= 0.1)
                            IsPreDry = false;
                        else
                            IsPreDry = true;
                    }                    
                    if (DryDays > 0)
                    {
                        double ProDW = DWCount * 1.0 / DryDays;
                        YearsCount += 1;
                        All += ProDW;
                    }
                }
                this.pData[month].PDW = All / YearsCount;
            }
        }
        /// <summary>
        /// 统计前一天为雨天当前天为雨天的概率
        /// </summary>
        private void CalProWW()
        {
            
            for (int month = 0; month < 12; month++)
            {
                int PreMonth = month - 1;
                double All = 0;
                int YearsCount =0;
                foreach (int year in WeatherData.Keys)
                {
                    int WetDays = 0;//这个月的雨天数
                    bool IsPreWet = false;
                    if (PreMonth < 0)
                    {
                        IsPreWet = true;
                    }
                    else if (WeatherData[year].Precipitation[PreMonth, this.daysInMonth[PreMonth] - 1] >= 0.1)  //从上个月的最后一天开始
                    {
                        IsPreWet = true;
                    }
                    int DDCount = 0;
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        if (IsPreWet)
                        {
                            WetDays += 1;
                            if (WeatherData[year].Precipitation[month, day] >= 0.1)
                            {
                                DDCount += 1;
                            }
                        }
                        if (WeatherData[year].Precipitation[month, day] >= 0.1)
                            IsPreWet = true;
                        else
                            IsPreWet = false;
                    }
                    if (WetDays > 0) //只有雨天不为0时，所有的才能算入
                    {
                        double ProWW = DDCount * 1.0 / WetDays;
                        All += ProWW;
                        YearsCount += 1;
                    }
                }
                this.pData[month].PWW = All / YearsCount;
            }
        }
        private void CalProWetDays()
        {
            int YearsCount = WeatherData.Keys.Count;
            for (int month = 0; month < 12; month++)
            {
                double All = 0;
                foreach (int year in WeatherData.Keys)
                {
                    int WetDays = 0;
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        if (WeatherData[year].Precipitation[month, day] > 0.1)
                            WetDays += 1;
                    }
                    All += WetDays / this.daysInMonth[month];
                }
                this.pData[month].PR_W = All / YearsCount;
            }
        }
        /// <summary>
        /// 计算按月统计值，计算以下15个统计值
        /// </summary>
        public void CalcMonthStat()
        {
            this.CalTempMax();
            this.CalTempMean();
            this.CalTempMeanSTD();
            this.CalTempMaxSTD();
            //this.CalSolarRadMean();
            this.CalWindSpeedMean();
            this.CalPreciMean();
            this.CalPreciSTD();
            this.CalDewPoint();
            this.CalMeanPCPMonSum();
            this.CalPCPDays();
            this.CalPCPSkew();
            this.CalProDW();
            this.CalProWW();
            this.CalProWetDays();            
        }

        //存放所有气象数据的结构，其中键为年份，值为一年中所有的值
        public SortedList<int, WeatherData> WeatherData = new SortedList<int, WeatherData>();
        //统放统计结果的结构
        public Statistics[] pData = new Statistics[12];
        /// <summary>
        /// 计算饱和水汽压
        /// </summary>
        /// <param name="airtemp"></param>
        /// <returns></returns>
        public double CalcSaturationVapP(double airtemp)
        {
            return Math.Exp((16.78*airtemp-116.9)/(airtemp+237.3));
        } 
        
        public void CalcPcpCoFactor(int pcpcode)
        {
            //long xrnd;
            //double rndm1, r6, rnm2, xlv, pcp, sum;

            //rndm1 = Random.NextDouble();

            //for (int mon = 0; mon < 12; mon++)
            //{
            //    if (pcpcode == 0)
            //    {
            //        r6 = 0;
            //        rnm2 = 0;
            //        xlv = 0;
            //        pcp = 0;
            //        sum = 0;
            //        r6 = MonData.pData[6, mon] / 6.0;
            //        for (int j = 0; j < 1000; j++)
            //        {
            //            rnm2 = Random.NextDouble();
            //            xlv = (GaussRandom() - r6) * r6 + 1;
            //            rndm1 = rnm2;
            //            xlv = (Math.Pow(xlv, 3) - 1.0) * 2 / MonData.pData[6, mon];
            //            pcp = xlv * MonData.pData[5, mon] + MonData.pData[15, mon];
            //            if (pcp < 0.01) pcp = 0.01;
            //            sum += pcp;
            //        }
            //        if (sum > 0)
            //            pcf[mon] = 1000 * MonData.pData[15, mon] / sum;
            //        else
            //            pcf[mon] = 1.0;
            //    }
            //}
        }
        
        /// <summary>
        /// 生成mon月strdate那天的降雨量
        /// </summary>
        /// <param name="mon"></param>
        /// <param name="strdate"></param>
        public void DailyPcp(int mon,DateTime strdate)
        {
            this.GenPcp_Skewed(strdate,true);
            throw new Exception("上面的程序还不是很完善");           
            
        }
        public void DailySlr(int mon, DateTime strdate)
        {
            //SlrGen(mon, strdate, true, 10, 30, 0.5, 0.85);
            throw new Exception("上面的程序还不是很完善");
            
        }
        public void DailyHmd(int mon, DateTime strdate)
        {
            this.HmdGen(mon, strdate);
            throw new Exception("上面的程序还不是很完善");
        }
        public void DailyWnd(int mon, DateTime strdate)
        {
            this.WndGen(mon, strdate);
            throw new Exception("上面的程序还不是很完善");
        }
        public void DailyTmp(int mon, DateTime strdate)
        {
            this.TempGen(mon, strdate, true, 10, 30, 0.5);
           throw new Exception("上面的程序还不是很完善");
        }
           
       
        public List<string> m_saIn, m_saOut;
        //CRandomv		randv;
        public double[] e = new double[3];
        /// <summary>
        /// 当前天的残差
        /// </summary>
        public double[] wgncur = new double[3];
        /// <summary>
        /// 前一天的残差
        /// </summary>
        public double[] wgnold = new double[3];
        private bool LoadPrecipitation = false;
        private bool LoadRelaHumidity = false;
        private bool LoadWindSpeed = false;
        private bool LoadMaxTemprature = false;
        private bool LoadMeanTemprature = false;
        private bool LoadDewPoint = false;
        private bool LoadSolarRadiation = false;

        /// <summary>
        /// 标准高斯分布的随机数
        /// </summary>
        /// <param name="rand1"></param>
        /// <param name="rand2"></param>
        /// <returns></returns>
        double GaussRandom()
        {
            double rand1 = Random.NextDouble();
            double rand2 = Random.NextDouble();
            return Math.Sqrt(-2.0 * Math.Log(rand1)) * Math.Cos(2 * Math.PI * rand2);
        }
        /// <summary>
        /// 读入降水数据
        /// </summary>
        public bool ReadRainfallData(string FileName)
        {
            StreamReader sr = new StreamReader(FileName);
            //先绕过两行
            sr.ReadLine();
            sr.ReadLine();
            while (!sr.EndOfStream)
            {
                string[] data= sr.ReadLine().Split(new char[]{' ','\t'},StringSplitOptions.RemoveEmptyEntries);
                int year=Convert.ToInt16( data[0].Substring(0, 4));
                if (this.WeatherData.ContainsKey(year)==false)
                {
                    this.WeatherData[year] = new WeatherData();
                }
                double[,] Data=this.WeatherData[year].Precipitation;
                int month = Convert.ToInt16(data[0].Substring(4,2))-1;
                for (int day = 0; day < this.daysInMonth[month]; day++)
                {
                    Data[month, day] = Convert.ToDouble(data[day + 1]);
                }
            }
            this.LoadPrecipitation = true;
            sr.Close();
            return true;
        }
        /// <summary>
        /// 读入最高温度数据
        /// </summary>
        /// <param name="FileName"></param>
        public bool ReadMaxTempData(string FileName)
        {
            StreamReader sr = new StreamReader(FileName);
            //先绕过两行
            sr.ReadLine();
            sr.ReadLine();
            while (!sr.EndOfStream)
            {
                string[] data = sr.ReadLine().Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
                int year = Convert.ToInt16(data[0].Substring(0, 4));
                if (this.WeatherData.ContainsKey(year) == false)
                {
                    this.WeatherData[year] = new WeatherData();
                }
                double[,] Data = this.WeatherData[year].MaxTemperature;
                int month = Convert.ToInt16(data[0].Substring(4, 2)) - 1;
                for (int day = 0; day < this.daysInMonth[month]; day++)
                {
                    Data[month, day] = Convert.ToDouble(data[day + 1])/10;
                }
            }
            sr.Close();
            this.LoadMaxTemprature = true;
            return true;
        }
        /// <summary>
        /// 读入平均温度数据
        /// </summary>
        /// <param name="FileName"></param>
        public bool ReadMeanTempData(string FileName)
        {
            StreamReader sr = new StreamReader(FileName);
            //先绕过两行
            sr.ReadLine();
            sr.ReadLine();
            while (!sr.EndOfStream)
            {
                string[] data = sr.ReadLine().Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
                int year = Convert.ToInt16(data[0].Substring(0, 4));
                if (this.WeatherData.ContainsKey(year) == false)
                {
                    this.WeatherData[year] = new WeatherData();
                }
                double[,] Data = this.WeatherData[year].MeanTemperature;
                int month = Convert.ToInt16(data[0].Substring(4, 2)) - 1;
                for (int day = 0; day < this.daysInMonth[month]; day++)
                {
                    Data[month, day] = Convert.ToDouble(data[day + 1]) / 10;
                }
            }
            sr.Close();
            this.LoadMeanTemprature = true;
            return true;
        }
        /// <summary>
        /// 读入日平均风速数据
        /// </summary>
        /// <param name="FileName"></param>
        public void ReadWindSpeed(string FileName)
        {
            StreamReader sr = new StreamReader(FileName);
            //先绕过两行
            sr.ReadLine();
            sr.ReadLine();
            while (!sr.EndOfStream)
            {
                string[] data = sr.ReadLine().Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
                int year = Convert.ToInt16(data[0].Substring(0, 4));
                if (this.WeatherData.ContainsKey(year) == false)
                {
                    this.WeatherData[year] = new WeatherData();
                }
                double[,] Data = this.WeatherData[year].WindSpeed;
                int month = Convert.ToInt16(data[0].Substring(4, 2)) - 1;
                for (int day = 0; day < this.daysInMonth[month]; day++)
                {
                    Data[month, day] = Convert.ToDouble(data[day + 1]);
                }
            }
            sr.Close();
            this.LoadWindSpeed = true;
            return;
        }
        /// <summary>
        /// 读入日平均相对湿度数据
        /// </summary>
        /// <param name="FileName"></param>
        public bool ReadRelaHumid(string FileName)
        {
            StreamReader sr = new StreamReader(FileName);
            //先绕过两行
            sr.ReadLine();
            sr.ReadLine();
            while (!sr.EndOfStream)
            {
                string[] data = sr.ReadLine().Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
                int year = Convert.ToInt16(data[0].Substring(0, 4));
                if (this.WeatherData.ContainsKey(year) == false)
                {
                    this.WeatherData[year] = new WeatherData();
                }
                double[,] Data = this.WeatherData[year].RelaHumidity;
                int month = Convert.ToInt16(data[0].Substring(4, 2)) - 1;
                for (int day = 0; day < this.daysInMonth[month]; day++)
                {
                    Data[month, day] = Convert.ToDouble(data[day + 1])/100;
                }
            }
            sr.Close();
            this.LoadRelaHumidity = true;
            return true;
        }
        
        /// <summary>
        /// 利用平均气温和相对湿度计算露点温度(根据SWAT理论推算)
        /// </summary>
        /// <param name="meanTemp">月平均气温</param>
        /// <param name="humidity">月平均相对湿度</param>
        /// <returns>露点温度</returns>
        public double CalcMeanDewPoint(double meanTemp,double humidity)
        {
            double k = Math.Log(humidity, Math.E) + (16.78 * meanTemp - 116.9) / (meanTemp + 237.3);
            return (237.3 * k + 116.9) / (16.78 - k);
        }
        /// <summary>
        /// 计算每日的平均露点温度(注意不是每月)
        /// </summary>
        public bool PrepareAllDewPoint()
        {
            if (this.LoadMeanTemprature == false || this.LoadRelaHumidity==false)
                throw new Exception("还没有导入所需要的数据");
            foreach (int year in this.WeatherData.Keys)
            {
                for (int month = 0; month < 12; month++)
                {
                    for (int day = 0; day < this.daysInMonth[month]; day++)
                    {
                        double meanTemp=this.WeatherData[year].MeanTemperature[month, day];
                        double humidity = this.WeatherData[year].RelaHumidity[month, day];
                        this.WeatherData[year].DewPoint[month, day] = this.CalcMeanDewPoint(meanTemp,humidity);
                    }
                }
            }
            this.LoadDewPoint = true;
            return true;
        }
        
    }
}