SimulationServer/Models/SimulationCommon/SeaSJ.cs

using OpenCvSharp;

namespace SimulationCommon;

public class WeatherResponse
{
    public Wind Wind { get; set; }
}

public class Wind
{
    public double Speed { get; set; }
    public double Deg { get; set; }
}

public class CurrentResponse
{
    public Current Current { get; set; }
}

public class Current
{
    public double Speed { get; set; }
    public double Direction { get; set; }
}

public class NCread
{
        
    //各数组来源文件'Text_readNC.cs'；时间范围：2024-06-04T00:00:00 ... 2024-06-04T23:00:00；
    //地点范围:落水点为中心100公里
    public float[] longitudeArray = new float[10];//经度一维数组来源自文件'Text_readNC.cs'
    public float[] latitudeArray = new float[7];//纬度一维数组来源自文件'Text_readNC.cs'
    public float[][][] u10Array ;//风的10米U（向东）分量三维数组来源自文件'Text_readNC.cs'
    public float[][][] v10Array;//风的10米V（向北）分量三维数组来源自文件'Text_readNC.cs'
    public float[][][] p140208Array;//海洋上空的自由对流速度三维数组来源自文件'Text_readNC.cs'
    public float[][][] mwdArray;//平均波向（单位：度；0度表示北方，90度表示东方）三维数组来源自文件'Text_readNC.cs'       
    public float[][][] hmaxArray;

}

public class SeaSJ
{
    // 生成符合正态分布的随机数
    public static double NormalDistributionRandom(double mu, double sigma)
    {
        Random rand = new Random();
        double u1 = 1.0 - rand.NextDouble(); // uniform(0,1] random doubles
        double u2 = 1.0 - rand.NextDouble();
        double normalRandom = Math.Sqrt(-2.0 * Math.Log(u1)) * Math.Sin(2.0 * Math.PI * u2); // random normal(0,1)
        return mu + sigma * normalRandom; // random normal(mean,stdDev^2)
    }

    // 交付类方法
    public static List<double[]> GetDrift(NCread nCread,double[] initialPosition, double dt, double totalTime)
    {
        // run 获取轨迹的函数
        // initialPosition --> 初始位置; dt --> 时间步长; totalTime --> 总时长
        List<double[]> trajectory = CalculateDriftTrajectory(nCread,initialPosition, dt, totalTime);
        return trajectory;
    }

    public static List<double[]> CalculateDriftTrajectory(NCread nCread,double[] initialPosition, double dt, double totalTime)
    {
        int timeSteps = (int)(totalTime / dt);
        List<double[]> trajectory = new List<double[]>();
        trajectory.Add(new double[] { initialPosition[0], initialPosition[1] });

        for (int t = 1; t < timeSteps; t++)
        {
            double[] currentPos = trajectory[t - 1];

            int outputTime = 3600 * t;
            // 动态获取当前位置的风力和洋流数据
            double[] windVelocity = GetWindVelocityFromAPI(nCread,currentPos[0], currentPos[1], outputTime);
            double[] currentVelocity = GetCurrentVelocityFromAPI(nCread,currentPos[0], currentPos[1], outputTime);

            // 计算漂移速度(m/s)
            double[] driftVelocity = {
                currentVelocity[0] + windVelocity[0],
                currentVelocity[1] + windVelocity[1]
            };

            // 更新位置（漂移速度单位是m/s，需要转换成经纬度的变化量）
            double[] newPosition = {
                currentPos[0] + (driftVelocity[0] * dt * 3600) / ((6371000 * Math.Cos(currentPos[0] * Math.PI / 180)) * 180 / Math.PI), // 纬度变化量
                currentPos[1] + (driftVelocity[1] * dt * 3600) / ((6371000 * Math.Cos(currentPos[1] * Math.PI / 180)) * 180 / Math.PI)  // 经度变化量
            };

            trajectory.Add(newPosition);
        }

        return trajectory;//（纬度，经度）
    }

    public static double[] GetWindVelocityFromAPI(NCread windNCread,double latitude, double longitude, double temptime)
    {
        float[] longitudeArray = windNCread.longitudeArray;//经度一维数组来源自文件'Text_readNC.cs'
        float[] latitudeArray = windNCread.latitudeArray;//纬度一维数组来源自文件'Text_readNC.cs'
        float[][][] u10Array = windNCread.u10Array;//风的10米U（向东）分量三维数组来源自文件'Text_readNC.cs'
        float[][][] v10Array = windNCread.v10Array;//风的10米V（向北）分量三维数组来源自文件'Text_readNC.cs'

        int longitudeNum = 0;
        int latitudeNum = 0;
        int temptimeNum = 0;
        //定义NC文件中读取不到的坐标点海洋数据信息（后续可采用插值法等）
        //经度连续化
        for (int i = 0; i < 9; i++)
        {
            if (longitude >= longitudeArray[i] && longitude < longitudeArray[i + 1])
            {
                longitude = longitudeArray[i];
                longitudeNum = i;
            }
        }

        //纬度连续化
        for (int i = 0; i < 6; i++)
        {
            if (latitude >= latitudeArray[i] && latitude < latitudeArray[i + 1])
            {
                latitude = latitudeArray[i];
                latitudeNum = i;
            }
        }
        
        //时间连续化
        for (int i = 0; i < 23 ; i ++)
        {
            if (temptime >= 3600 * i && temptime < (3600 * (i + 1)))
            {
                temptimeNum = i;
            }
        }

        double windX = (double)u10Array[temptimeNum][latitudeNum][longitudeNum];
        double windY = (double)v10Array[temptimeNum][latitudeNum][ longitudeNum];

        return new double[] { windX, windY };
    }

    public static double[] GetCurrentVelocityFromAPI(NCread CurrentNCread,double latitude, double longitude, double temptime)
    {
        float[] longitudeArray = CurrentNCread.longitudeArray;//经度一维数组来源自文件'Text_readNC.cs'
        float[] latitudeArray = CurrentNCread.latitudeArray;//纬度一维数组来源自文件'Text_readNC.cs'
        float[][][] u10Array = CurrentNCread.u10Array;//风的10米U（向东）分量三维数组来源自文件'Text_readNC.cs'
        float[][][] v10Array = CurrentNCread.v10Array;//风的10米V（向北）分量三维数组来源自文件'Text_readNC.cs'
        float[][][] p140208Array = CurrentNCread.p140208Array;//海洋上空的自由对流速度三维数组来源自文件'Text_readNC.cs'
        float[][][] mwdArray = CurrentNCread.mwdArray;//平均波向（单位：度；0度表示北方，90度表示东方）三维数组来源自文件'Text_readNC.cs'

        int longitudeNum = 0;
        int latitudeNum = 0;
        int temptimeNum = 0;
        //定义NC文件中读取不到的坐标点海洋数据信息（后续可采用插值法等）
        //经度连续化
        for (int i = 0; i < 9; i++)
        {
            if (longitude >= longitudeArray[i] && longitude < longitudeArray[i + 1])
            {
                longitude = longitudeArray[i];
                longitudeNum = i;
            }
        }

        //纬度连续化
        for (int i = 0; i < 6; i++)
        {
            if (latitude >= latitudeArray[i] && latitude < latitudeArray[i + 1])
            {
                latitude = latitudeArray[i];
                latitudeNum = i;
            }
        }
        
        //时间连续化
        for (int i = 0; i < 23; i++)
        {
            if (temptime >= 3600 * 1 && temptime < 3600 * (i + 1))
            {
                temptimeNum = i;
            }
        }
        double currentSpeed = (double)p140208Array[temptimeNum][latitudeNum][ longitudeNum];
        double currentDirection = (double)mwdArray[temptimeNum][latitudeNum][ longitudeNum];

        double currentDirectionInRadians = currentDirection * (Math.PI / 180);
        double currentX = currentSpeed * Math.Cos(currentDirectionInRadians);
        double currentY = currentSpeed * Math.Sin(currentDirectionInRadians);

        return new double[] { currentX, currentY };
    } 
    
    //海浪高度获取
   public static double GetWaveHeightFromAPI(NCread WaveNCread,double latitude, double longitude, double temptime)//temptime:仿真时间（秒）
    {
        float[] longitudeArray = WaveNCread.longitudeArray;//经度一维数组来源自文件'Text_readNC.cs'
        float[] latitudeArray = WaveNCread.latitudeArray;//纬度一维数组来源自文件'Text_readNC.cs'
        float[][][] hmaxArray = WaveNCread.hmaxArray;//风的10米U（向东）分量三维数组来源自文件'Text_readNC.cs'

        int longitudeNum = 0;
        int latitudeNum = 0;
        int temptimeNum = 0;

        //定义NC文件中读取不到的坐标点海洋数据信息（后续可采用插值法等）
        //经度连续化
        for (int i = 0; i < 9; i++)
        {
            if (longitude >= longitudeArray[i] && longitude < longitudeArray[i + 1])
            {
                longitude = longitudeArray[i];
                longitudeNum = i;
            }
        }

        //纬度连续化
        for (int i = 0; i < 6; i++)
        {
            if (latitude >= latitudeArray[i] && latitude < latitudeArray[i + 1])
            {
                latitude = latitudeArray[i];
                latitudeNum = i;
            }
        }

        //时间连续化
        for (int i = 0; i < 23 ; i ++)
        {
            if (temptime >= 3600 * i && temptime < (3600 * (i + 1)))
            {
                temptimeNum = i;
            }
        }

        double WaveHeight = (double)hmaxArray[temptimeNum][ latitudeNum][ longitudeNum];

        return WaveHeight;
    }
    
    // 交付部分方法代码
    public static List<double[]> getminEnclosingRect(List<double[]> latLonList)
    {
        // 转换经纬度为墨卡托投影坐标，并添加到点集合
        List<Point2f> pointList = new List<Point2f>();
        foreach (var latLon in latLonList)
        {
            double x = Rectangular_Area_Search_Function.MokatuoLat(latLon[1]);
            double y = Rectangular_Area_Search_Function.MokatuoLon(latLon[0]);
            pointList.Add(new Point2f((float)y, (float)x));
        }

        // 获取凸包
        Point2f[] convexHull = Rectangular_Area_Search_Function.GetConvexHull(pointList);
        
        //获取凸包各点经纬度坐标
        List<double[]> hullPoint = new List<double[]>();
        for (int num = 0; num < convexHull.Length - 1; num++)
        {
            hullPoint.Add(new double[] { convexHull[num].Y, convexHull[num].X });
        }

        List<double[]> hullPointLatLon = new List<double[]>();
        foreach (var point in hullPoint)
        {
            double pointLat = Rectangular_Area_Search_Function.RMokatuoLat(point[1]);
            double pointLon = Rectangular_Area_Search_Function.RMokatuoLon(point[0]);
            hullPointLatLon.Add(new double[] { pointLat, pointLon });
        }
        
        // 计算最小包围矩形
        List<double[]> minEnclosingRect = Rectangular_Area_Search_Function.MinEnclosingRectangle(convexHull);

        // 最小包围矩形顶点经纬度坐标
        List<double[]> startPoint = new List<double[]>();
        foreach (var minEnclosingRectPoint in minEnclosingRect)
        {
            double lat = Rectangular_Area_Search_Function.RMokatuoLat(minEnclosingRectPoint[0]);
            double lon = Rectangular_Area_Search_Function.RMokatuoLon(minEnclosingRectPoint[1]);
            startPoint.Add(new double[] { lat, lon });
        }

        return startPoint;
    }
}