【Matlab 根据纬度求一天光周期时长 daylength】

daylength is a function of latitude and day of year

因为工作中用到了计算一天光周期，所以参考网上零零散散的内容，用matlab写了一下：

%     y=year
%     m=month
%     d=day
%     formatOut = 'yyyy-mm-dd';
%     datas=datestr(datenum(y,m,d),formatOut,'local');
%     yDay=sum(eomday(y, 1:m))+d;

function DL= LatToDayLength(lat,year,month,day)

    %get lat into radians
    gamma=lat/180.0*pi;
    %convert data into an angle
    yDay=sum(eomday(year,1:month))+day;
    psi=yDay/365.0*2*pi;
   % calculate solar declination
    solarDec=(0.39637-22.9133*cos(psi)+4.02543*sin(psi)-0.38720*cos(2*psi)+0.05200*sin(2*psi))*pi/180.0;
    r=-tan(gamma)*tan(solarDec);
if r<= -1
    DL=24.0;
else
    if abs(r)<1
        DL=24*acos(r)/pi;
    else
        DL=0;
    end
end