I annotated the gravitational force and it worked, but I got an error when I added it.The error message is org.orekit.errors.OrekitException: invalid parameter eccentricity: -0 not in range [0, ∞].Here are some code snippets I wrote. Could you please help me check if there are any issues?
public void SatelliteVisibility1() throws Exception {
//加载固定数据
Dataloading dataloading = new Dataloading();
dataloading.dataLoade();
//时间设置
String startTime = “2024-12-23 10:00:00”;
String endTime = “2024-12-23 15:00:00”;
SimpleDateFormat format = new SimpleDateFormat(“yyyy-MM-dd HH:mm:ss”);
Date StartTime = format.parse(startTime);
Date EndTime = format.parse(endTime);
// Initial date in UTC time scale:定义时间类型
final TimeScale utc = TimeScalesFactory.getUTC();
//转换开始/结束时间类型
AbsoluteDate start = new AbsoluteDate(StartTime, utc);
AbsoluteDate end = new AbsoluteDate(EndTime, utc);
//转换轨道历元时间类型
final AbsoluteDate initialDate = new AbsoluteDate(2024, 12, 23, 02, 38, 36.000, utc);
//轨道六根数
final double a = 6893678.08724067; // 半长轴 (m)
final double e = 0.00153912; // 离心率
final double i = FastMath.toRadians(97.60317163); // 轨道倾角 (弧度)
final double omega = FastMath.toRadians(65.73638295); // 近地点幅角 (弧度)
final double raan = FastMath.toRadians(69.24863943); // 升交点赤经 (弧度)
final double meanAnomaly = FastMath.toRadians(55.36150646); // 平近点角 (弧度)
// Inertial frame
final Frame inertialFrame = FramesFactory.getEME2000();
// gravitation coefficient:定义地球引力系数
final double mu = 3.986004415e+14;
// Orbit construction as Keplerian:定义卫星轨道类型
final Orbit initialOrbit = new KeplerianOrbit(a, e, i, omega, raan, meanAnomaly, PositionAngleType.MEAN,inertialFrame, initialDate, mu);
// Initialize state:定义卫星的初始状态
final SpacecraftState initialState = new SpacecraftState(initialOrbit, 31);
/–定义propagator并插入到卫星上–/
//使用积分器是 8(5,3) 阶的嵌入式 Runge-Kutta 积分器,用于局部外推模式
// steps limits
final double minStep = 0.001;
final double maxStep = 1000;
final double initStep = 60;
// error control parameters (absolute and relative)
final double positionError = 10.0;
final double tolerances = NumericalPropagator.tolerances(positionError, initialOrbit, initialOrbit.getType());
// set up mathematical integrator
AdaptiveStepsizeIntegrator integrator =
new DormandPrince853Integrator(minStep, maxStep, tolerances[0], tolerances[1]);
integrator.setInitialStepSize(initStep);
//使用龙格库塔积分器
//final double stepSize = 10;
//final AbstractIntegrator integrator = new ClassicalRungeKuttaIntegrator(stepSize);
final NumericalPropagator numericalPropagator = new NumericalPropagator(integrator);
// Earth and frame:定义地球框架
final Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
final BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
Constants.WGS84_EARTH_FLATTENING, earthFrame);
/–力学模型设置–/
//设置重力场模型
int degree = 21;
int order = 21;
NormalizedSphericalHarmonicsProvider gravProvider = GravityFieldFactory.getNormalizedProvider(degree, order); // WGS84 potential coefficients
ForceModel gravModel = new HolmesFeatherstoneAttractionModel(earthFrame, gravProvider);
numericalPropagator.addForceModel(gravModel);
//设置潮汐力模型
// Solid Earth Tides
SolidTides tides = new SolidTides(earthFrame,
Constants.EIGEN5C_EARTH_EQUATORIAL_RADIUS, Constants.EIGEN5C_EARTH_MU,
TideSystem.ZERO_TIDE,
IERSConventions.IERS_2010,
TimeScalesFactory.getUT1(IERSConventions.IERS_2010, false),
CelestialBodyFactory.getSun(), CelestialBodyFactory.getMoon());
numericalPropagator.addForceModel(tides);
// Relativity
Relativity rel = new Relativity(Constants.EIGEN5C_EARTH_MU);
numericalPropagator.addForceModel(rel);
//设置第三方引力模型
numericalPropagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun()));
numericalPropagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon()));
// // Solar Radiation Pressure Model
// if (!Double.isNaN(coefAbsorb) && !Double.isNaN(coefReflect)) {
// RadiationSensitive satSRPConfig = new IsotropicRadiationClassicalConvention(crossSection, coefAbsorb,
// coefReflect);
// SolarRadiationPressure srp = new SolarRadiationPressure(CelestialBodyFactory.getSun(),
// Constants.EIGEN5C_EARTH_EQUATORIAL_RADIUS, satSRPConfig);
// srp.addOccultingBody(CelestialBodyFactory.getMoon(),
// Constants.MOON_EQUATORIAL_RADIUS);
// numericalPropagator.addForceModel(srp);
// }
//Drag Model:设置大气阻力模型
CssiSpaceWeatherData inputParameters = new CssiSpaceWeatherData(“SpaceWeather-All-v1.2.txt”);
Atmosphere atmosphere = new NRLMSISE00(inputParameters, CelestialBodyFactory.getSun(), earth);
ForceModel dragModel = new DragForce(atmosphere, new IsotropicDrag(0.62396452, 2.2));//迎风面与阻力系数
numericalPropagator.addForceModel(dragModel);
// Initialize propagation:初始化传播器
numericalPropagator.setInitialState(initialState);
numericalPropagator.setOrbitType(initialOrbit.getType());
numericalPropagator.setPositionAngleType(PositionAngleType.MEAN);
//定义地面站:名称、经度、纬度、高度(m)、起跟仰角
final String facilityName = “JH0751”;
final double lon = FastMath.toRadians(82.865);
final double lat = FastMath.toRadians(44.560);
final double altitude = 360;
final double minElevation = FastMath.toRadians(3.0);
final GeodeticPoint facility = new GeodeticPoint(lat, lon, altitude);
final TopocentricFrame sta1Frame = new TopocentricFrame(earth, facility, facilityName);
//定义事件处理器,用于检测可见弧段时间
final double maxcheck = 60.0;
final double threshold = 0.01;
final EventDetector sta1Visi =
new ElevationDetector(maxcheck, threshold, sta1Frame)
.withConstantElevation(minElevation)
.withHandler((s, detector, increasing) → {
if (increasing) {
System.out.println(“地面站名称: “+((ElevationDetector) detector).getTopocentricFrame().getName()
+” 跟踪开始时间: " + s.getDate()
+” 入境方位角: “+FastMath.toDegrees(((ElevationDetector) detector).getTopocentricFrame().getAzimuth(s.getPosition(), s.getFrame(), s.getDate()))
+” 入境俯仰角: "+FastMath.toDegrees(((ElevationDetector) detector).getTopocentricFrame().getElevation(s.getPosition(), s.getFrame(), s.getDate()))
);
} else {
System.out.println(“地面站名称: “+((ElevationDetector) detector).getTopocentricFrame().getName()
+” 跟踪结束时间: " + s.getDate()
+” 出境方位角: “+FastMath.toDegrees(((ElevationDetector) detector).getTopocentricFrame().getAzimuth(s.getPosition(), s.getFrame(), s.getDate()))
+” 出境俯仰角: "+FastMath.toDegrees(((ElevationDetector) detector).getTopocentricFrame().getElevation(s.getPosition(), s.getFrame(), s.getDate()))
);
}
return Action.CONTINUE; // Continue processing the event
});
// Add event to be detected
numericalPropagator.addEventDetector(sta1Visi);
//定义事件处理器,用于检测最大仰角
ElevationExtremumDetector extremumDetector = new ElevationExtremumDetector(sta1Frame);
final EventSlopeFilter maxElevationDetector = new EventSlopeFilter<>(extremumDetector,
FilterType.TRIGGER_ONLY_DECREASING_EVENTS);
EventHandler handler = (s, detector, increasing) → {
double elevation = extremumDetector.getElevation(s);
//将仰角为负的点过滤掉
if (elevation<FastMath.toRadians(3.0)) {
return Action.CONTINUE;
}
//判断卫星速度矢量方向向北分量是否为正
//if (s.getPVCoordinates().getVelocity().getNorm()>sta1Frame.getNorth().getNorm())
if(sta1Frame.getNorth().dotProduct(s.getPVCoordinates().getVelocity())<0)
{
System.out.println(" 降轨圈次过顶时间:"+s.getDate()+" 最大仰角:"+FastMath.toDegrees(elevation));
}
else{
System.out.println(" 升轨圈次过顶时间:"+s.getDate()+" 最大仰角:"+FastMath.toDegrees(elevation));
}
return Action.CONTINUE;
};
numericalPropagator.addEventDetector(maxElevationDetector.withHandler(handler));
// Propagate from the initial date to the first raising or for the fixed duration
final SpacecraftState state = numericalPropagator.propagate(start,end);
}