测站预报（基于OREKIT）

public static ArrayList<Double> stationPropagateAngle(final String Line1, final String Line2, final String line1, final String line2, final TimeScale utc,
                                                          final double longitude, final double latitude, final double altitude,
                                                          final double maxcheck, final double threshold, final double elevation,
                                                          final AbsoluteDate initialDate,  final double duration,double step ) {
        ArrayList<Double> angle = new ArrayList<>();
        ArrayList<Vector3D> t1 = new ArrayList<>();
        ArrayList<Vector3D> t2 = new ArrayList<>();
        ArrayList<Vector3D> st = new ArrayList<>();
        //常量
        final double mu = 3.986004415e+14; // gravitation coefficient 重力常数
        final Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
        final BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
                Constants.WGS84_EARTH_FLATTENING,
                earthFrame);
        final Frame inertialFrame = FramesFactory.getEME2000(); // inertial frame for orbit definition 用于轨道定义的惯性系

        //积分器
        final TLE tle1 = new TLE(Line1, Line2);
        final Propagator propagator = TLEPropagator.selectExtrapolator(tle1);
        final Propagator propagator1 = TLEPropagator.selectExtrapolator(tle1);//主星
        final TLE tle2 = new TLE(line1, line2);
        final Propagator propagator2 = TLEPropagator.selectExtrapolator(tle2);

        //测站
        final GeodeticPoint station = new GeodeticPoint(latitude, longitude, altitude);
        final TopocentricFrame sta1Frame = new TopocentricFrame(earth, station, "station1");

        final AbsoluteDate[] dates = {new AbsoluteDate(), new AbsoluteDate()};

        //检测事件
        final EventDetector sta1Visi =
                new ElevationDetector(maxcheck, threshold, sta1Frame).
                        withConstantElevation(elevation).
                        withHandler((s, detector, increasing) -> {
                            //AbsoluteDate BeginTime =
                            if (increasing) {
                                dates[0] = s.getDate();
                            } else {
                                dates[1] = s.getDate();
                            }
                            return increasing ? Action.CONTINUE : Action.STOP;
                        });

        //添加要检测的事件
        propagator.addEventDetector(sta1Visi);
        final SpacecraftState finalState = propagator.propagate(initialDate.shiftedBy(duration));

        //计算角度
        AbsoluteDate time = dates[0];
        do {
            //测站坐标系转J2000
            Transform statoJ2000 = sta1Frame.getTransformTo(FramesFactory.getEME2000(), time);
            Transform tle1oJ2000 = propagator1.getFrame().getTransformTo(FramesFactory.getEME2000(),time);
            Transform tle2oJ2000 = propagator2.getFrame().getTransformTo(FramesFactory.getEME2000(),time);
            Vector3D sta = statoJ2000.transformPosition(new Vector3D(latitude, longitude, altitude));
            Vector3D p1 = tle1oJ2000.transformPosition(propagator1.propagate(time).getPVCoordinates().getPosition());
            Vector3D p2 = tle2oJ2000.transformPosition(propagator2.propagate(time).getPVCoordinates().getPosition());
            Vector3D v1 = p1.subtract(sta);
            Vector3D v2 = p2.subtract(sta);
            double p = v1.dotProduct(v2) / (v1.getNorm() * v2.getNorm());
            double Angle = Math.toDegrees(Math.acos(p));
            angle.add(Angle);
            t1.add(p1);
            t2.add(p2);
            st.add(sta);
            time = time.shiftedBy(step);
        } while (time.isBetweenOrEqualTo(dates[0],dates[1]));
        return angle;
    }