日潜在蒸散发彭曼公式Matlab

%
%输入：
%   tasmax：日最高气温
%   tasmin：日最低气温
%   rhs：近地相对湿度
%   wind0：近地风速
%   rsds：短波辐射
%   lat：纬度
%   alt：海拔
%   time：时间
%
function PE = CalculatePE（tasmax tasmin rhs wind rsds lat alt time）
    % inputs and units changing
    %Tav=tas-273.16;  % tas （surface air tempK） （1K = （°C） + 273.16）
    Tmax=tasmax-273.16; % tasmax
    Tmin=tasmin-273.16; % tasmin
    rh=rhs/100;      % rhs （near-surface relative humidity %）
    %prec=pr*86400;  % pr （precipitation kg/m2/s） （*86400=mm/day）
    wind0=wind;        %  wind（near-sruface wind speed at 10mm/s） 
    Rs=rsds * 0.0864;    % rsds （shortwave downwelling w/m2） （0.0864convert from W/m2 to MJ/m2/day）;
    %Rlds=rlds * 0.0864;  % rlds （longwave downwelling w/m2） （0.0864convert from W/m2 to MJ/m2/day）;
    %alt=0;  %altitude （m）
    %lat=0;  %latitude （degree）
    %DOY=1;  %  DOY 1 to 365（366）
    %time=2013-3-2;  % time convert to DOY  by the followed:
    % time to DOY: matlab time vector （midpoint） （HAS TO BE A COLUMN VECTOR!）
    %dt = median（diff（time））;
    %time_vec = datevec（time）;   

    %tmean 
    T=（Tmax+Tmin）/2.0;       
    
    %esea  饱和水汽压和实际水汽压的计算
    etx=0.6108*exp（（17.27*Tmax）./（Tmax+237.3））;
    etn=0.6108*exp（（17.27*Tmin）./（Tmin+237.3））;
    es=（etx+etn）/2.0;
    ea=es.*rh;    
    
    %slope饱和水汽压曲线斜率计算
    slope=（4098*0.6108*exp（（17.27*T）./（T+237.3）））./（（T+237.3）.^2）;
    
    %Rn净辐射计算