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

namespace GeoFly
{
    public class SolarRadiation
    {
        /// <summary>
        /// 地球自转角速度，其值为0.2618弧度/h
        /// </summary>
        double AV = 0.2618;
        /// <summary>
        /// 太阳常数(MJ/M^2/J)
        /// </summary>
        double ISC = 4.921;
        public SolarRadiation(int year, int dn, double lat, double shour)
        {
        }
        /// <summary>
        /// 
        /// </summary>
        /// <param name="year"></param>
        /// <param name="dn"></param>
        /// <param name="lat"></param>
        public SolarRadiation(int year, int dn, double lat)
        {
        }
        /// <summary>
        /// 根据日期太阳直射在地球上的纬度(返回弧度值)
        /// </summary>
        /// /// <param name="dn">Julie历日期</param>
        /// <returns></returns>
        public double SolarDeclination(int dn)
        {
            double value = 0.4 * Math.Sin(2 * Math.PI * (dn - 82) / 365);
            return Math.Asin(value);	
        }
        /// <summary>
        /// 计算eccentricity correction factor
        /// </summary>
        /// <param name="dn">Julie历日期</param>
        /// <returns></returns>
        public double EccentricityCF(int dn)
        {            
	        return 1+0.033*Math.Cos(2*Math.PI*dn/365.0);
        }
        /// <summary>
        /// the time sun is up,the time interval is between sunrise to solar noon
        /// 
        /// </summary>
        /// <returns></returns>
        public double SunRiseTime(int dn,double Latitude)
        {
            return Math.Acos(-Math.Tan(this.SolarDeclination(dn)) *Math.Tan(Latitude*Math.PI/180.0)) / this.AV;
        }
        /// <summary>
        /// the time sun is down,the time interval is between solar noon to sunset
        /// </summary>
        /// <returns></returns>
        public double SunSetTime(int dn,double Latitude)
        {
            return -Math.Acos(-Math.Tan(this.SolarDeclination(dn)) * Math.Tan(Latitude * Math.PI / 180.0)) / this.AV;
        }
        /// <summary>
        /// the total day length at latitudes between 66.5 and -66.5
        /// </summary>
        /// <returns></returns>
        public double DayLength(int dn,double Latitude)
        {
            return 2*Math.Acos(-Math.Tan(this.SolarDeclination(dn)) * Math.Tan(Latitude * Math.PI / 180.0)) / this.AV;
        }
        /// <summary>
        /// the height between the sun and a horizontal surface on the earth's surface
        /// </summary>
        /// <returns></returns>
        public double SunHeight()
        {
            //double dret;
            //double ha = AV * dSHour;		//the hour angle
            //dret =Math.Sin(dSd) *Math.Sin(dLat*Math.PI/180.0) + Math.Cos(dSd) *Math.Cos(dLat*Math.PI/180.0) * Math.Cos(ha);
            //return dret;
            throw new Exception("asdfsafsdf");
        }
        /// <summary>
        /// 在给定纬度上给定日期的地球表面上的日天文辐射量(MJ/m^2/day)
        /// </summary>
        /// <returns></returns>
        public double ExtraTerrRad(int dn, double Latitude)
        {
            double tsr = SunRiseTime(dn, Latitude);
            double E0 = this.EccentricityCF(dn);
            double delta=this.SolarDeclination(dn);
            double phi = Latitude * Math.PI / 180;
            double value=this.AV*this.SunRiseTime(dn,Latitude);
            double dd = value * Math.Sin(delta) * Math.Sin(phi) + Math.Cos(delta) * Math.Cos(phi) * Math.Sin(value);
            return (24.0 / Math.PI) * ISC * E0 * dd;
           
        }
        /// <summary>
        ///晴朗天空下的太阳辐射，即最大可能的日太阳辐射
        /// </summary>
        /// <param name="dn">Julie历日期</param>
        /// <param name="Latitude">当地纬度(角度值)</param>
        /// <param name="dscale">晴朗大气对辐射的减弱因子</param>
        /// <returns></returns>
        public double RealSolarRadMax(int dn,double Latitude, double dscale)
        {            
            return this.ExtraTerrRad(dn, Latitude) * dscale;
        }

        /// <summary>
        ///  the real solar radiation under real cloud conditions
        /// may be the extraterrestrial radiation or the ideal atmospherical radiation
        ///		or the clear sky total solar radiation;
        /// </summary>
        /// <param name="slrpercent">sunlight percentage(日照百分率)</param>
        /// <param name="a">the experiential coefficient ie: y = a + b * s，默认为0.16</param>
        /// <param name="b">the experiential coefficient ie: y = a + b * s，默认为0.65</param>
        /// <param name="dabsort">默认为0.85</param>
        /// <returns></returns>
        public double RealSolarRad(int dn,double Latitude, double slrpercent, double a, double b, double dabsorb)
        {
            double Q0 = RealSolarRadMax(dn, Latitude, dabsorb);
            return Q0 * (a + b * slrpercent);
        }
        /// <summary>
        /// 同上，只是默认a=0.16,b=0.65,dabsort=0.85
        /// </summary>
        /// <param name="slrpercent"></param>
        /// <returns></returns>
        public double RealSolarRad(int dn,double Latitude,double slrpercent)
        {
            return RealSolarRad(dn,Latitude,slrpercent, 0.16, 0.65, 0.85);
        }
        /// <summary>
        /// 判断输入年份是否是闰年
        /// </summary>
        /// <param name="year">需判断年份</param>
        /// <returns>true 是闰年；false不是闰年</returns>
        public bool IsLeapYear(int year)
        {
            if (year % 4 == 0)			//Is it a century?
            {
                if (year % 100 == 0)    // If a century, must be evenly divisible by 400.
                {
                    if (year % 400 == 0)
                        return true;
                    else
                        return false;
                }
                else
                    return true;
            }
            else
                return false;
        }
    }
}