package com.dy.pipIrrModel.modelCalculate;
|
|
import net.objecthunter.exp4j.Expression;
|
import net.objecthunter.exp4j.ExpressionBuilder;
|
|
/**
|
* @Author: liurunyu
|
* @Date: 2025/8/18 14:51
|
* @Description 哈格里夫斯（Hargreaves）公式 计算作物实际蒸散量（mm/day）
|
*/
|
public class Hargreaves {
|
|
/**
|
* 计算弧度
|
* 弧度=角度× π/180
|
* @param lat 地理纬度
|
* @return
|
*/
|
public static Double rad(Double lat) {
|
Expression expression = new ExpressionBuilder("x * pi / 180")
|
.variables("x")
|
.build();
|
// 设置变量值
|
expression.setVariable("x", lat);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
/**
|
* 太阳磁偏角
|
* @param dayIndex 为年内某天的日序数（比如1月1日为1，12月31日为365或366）
|
* @return
|
*/
|
public static Double sunMagnetismAngular(Integer dayIndex){
|
Expression expression = new ExpressionBuilder("0.409 * sin(((2 * pi) / 365) * x - 1.39)")
|
.variables("x")
|
.build();
|
// 设置变量值
|
expression.setVariable("x", dayIndex);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
/**
|
* 日地间相对距离的倒数
|
* @param dayIndex 为年内某天的日序数（比如1月1日为1，12月31日为365或366）
|
* @return
|
*/
|
public static Double sunEarthDistance(Integer dayIndex){
|
Expression expression = new ExpressionBuilder("1 + (0.033 * cos(((2 * pi) / 365) * x))")
|
.variables("x")
|
.build();
|
// 设置变量值
|
expression.setVariable("x", dayIndex);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
|
/**
|
* 太阳时角
|
* 0.409 * sin(2 * pi / 365 * x - 1.39)
|
* @param fai 地理弧度
|
* @param sunMagnetismAngular 太阳磁偏角
|
* @return
|
*/
|
public static Double sunTimeAngular(Double fai, Double sunMagnetismAngular){
|
Expression expression = new ExpressionBuilder("acos(-tan(x) * tan(y))")
|
.variables("x", "y")
|
.build();
|
// 设置变量值
|
expression.setVariable("x", fai);
|
expression.setVariable("y", sunMagnetismAngular);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
/**
|
* 计算天顶辐射
|
* @param sunEarthDistance 日地间相对距离的倒数
|
* @param sunTimeAngular 太阳时偏角
|
* @param fai 地理弧度
|
* @param sunMagnetismAngular 太阳磁偏角
|
* @return
|
*/
|
public static Double zenithRadiation(Double sunEarthDistance,
|
Double sunTimeAngular,
|
Double fai,
|
Double sunMagnetismAngular){
|
Expression expression = new ExpressionBuilder("((24 * 60) / pi) * 0.082 * a * ((b * sin(c) * sin(d)) + (cos(c) * cos(d) * sin(b)))")
|
.variables("a", "b", "c", "d")
|
.build();
|
// 设置变量值
|
expression.setVariable("a", sunEarthDistance);
|
expression.setVariable("b", sunTimeAngular);
|
expression.setVariable("c", fai);
|
expression.setVariable("d", sunMagnetismAngular);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
/**
|
* 格里夫斯（Hargreaves）公式计算作物蒸散量（mm/day）
|
* @param kc 作物系数
|
* @param maxT 一日内最高温度
|
* @param maxT 一日内最高温度
|
* @param zenithRadiation 天顶辐射
|
* @return
|
*/
|
public static Double ET0(Double kc, Double maxT, Double minT, Double zenithRadiation){
|
Expression expression = new ExpressionBuilder("x * (0.0023 * (((a + b) / 2) + 17.8) * ((a + b)^(1/2)) * c * 0.408)")
|
.variables("x", "a", "b", "c")
|
.build();
|
// 设置变量值
|
expression.setVariable("x", kc);
|
expression.setVariable("a", maxT);
|
expression.setVariable("b", minT);
|
expression.setVariable("c", zenithRadiation);
|
Double result = expression.evaluate();
|
return result ;
|
}
|
|
public static void main(String[] args) {
|
Double lat = 38.561976140977116 ;//纬度
|
Integer dayIndex = 180 ;//年内日序数（比如1月1日为1，12月31日为365或366）
|
Double kc = 0.41 ;//作物系数
|
Double maxT = 40.1 ;//一日内最高温度
|
Double minT = 28.1 ;//一日内最低温度
|
|
Double fai = rad(lat);
|
System.out.println(fai);
|
|
Double sunMagnetismAngular = sunMagnetismAngular(dayIndex);
|
System.out.println(sunMagnetismAngular);
|
|
Double sunEarthDistance = sunEarthDistance(dayIndex);
|
System.out.println(sunEarthDistance);
|
|
Double sunTimeAngular = sunTimeAngular(fai, sunMagnetismAngular);
|
System.out.println(sunTimeAngular);
|
|
Double zenithRadiation = zenithRadiation(sunEarthDistance, sunTimeAngular, fai, sunMagnetismAngular);
|
System.out.println(zenithRadiation);
|
|
Double et0 = ET0(kc, maxT, minT, zenithRadiation);
|
System.out.println(et0);
|
}
|
|
}
|
