C语言版本路径压缩算法

line_vacuate/studio_geo_algo_c.c

@@ -0,0 +1,168 @@
+/**
+  ******************************************************************************
+  * @file           : studio_geo_algo_c.c
+  * @author         : wangyingjie
+  * @brief          : None
+  * @attention      : None
+  * @date           : 2025/5/11
+  ******************************************************************************
+  */
+
+#include "studio_geo_algo_c.h"
+
+
+double point_line_dist_square_c(const studio_point_c *p, const studio_point_c *start, const studio_point_c *end)
+{
+    double A = p->x - start->x;
+    double B = p->y - start->y;
+    double C = end->x - start->x;
+    double D = end->y - start->y;
+
+    double dot = A * C + B * D;
+    double len_sq = C * C + D * D;
+    double param = len_sq == 0 ? -1 : dot / len_sq;
+
+    double xx, yy;
+
+    if (param < 0)
+    {
+        xx = start->x;
+        yy = start->y;
+    }
+    else if (param > 1)
+    {
+        xx = end->x;
+        yy = end->y;
+    }
+    else
+    {
+        xx = start->x + param * C;
+        yy = start->y + param * D;
+    }
+
+    double dist_square = (xx - p->x) * (xx - p->x) + (yy - p->y) * (yy - p->y);
+    return dist_square;
+}
+
+// Douglas-Peucker算法，用于简化多边形或曲线
+void douglas_peucker_c(const studio_line_c *points, size_t start, size_t end, double epsilon, unsigned int *indices, unsigned int *index_count)
+{
+    // 计算最大距离
+    double max_dist = 0;
+    size_t index = start;
+    double epsilon_square = epsilon * epsilon;
+
+    // 遍历从start+1到end-1的点
+    for (size_t i = start + 1; i < end; ++i)
+    {
+        // 计算点到线段的距离的平方
+        double dist_square = point_line_dist_square_c(&points->data[i], &points->data[start], &points->data[end]);
+        // 如果距离大于最大距离，则更新最大距离和索引
+        if (dist_square > max_dist)
+        {
+            index = i;
+            max_dist = dist_square;
+        }
+    }
+
+    // 如果最大距离大于epsilon的平方，并且索引不是start和end，则递归调用Douglas-Peucker算法
+    if (max_dist > epsilon_square && index != start && index != end)
+    {
+        douglas_peucker_c(points, start, index, epsilon, indices, index_count);
+        douglas_peucker_c(points, index, end, epsilon, indices, index_count);
+    }
+    // 否则，将start和end的索引添加到indices数组中
+    else
+    {
+        indices[*index_count] = start;
+        (*index_count)++;
+        if (start != end)
+        {
+            indices[*index_count] = end;
+            (*index_count)++;
+        }
+    }
+}
+
+void remove_duplicates(unsigned int **indices, unsigned int *size) {
+    if (*size == 0) return;
+    // 临时数组用于存储去重后的数据
+    unsigned int *unique_indices = (unsigned int *)malloc(*size * sizeof(unsigned int));
+    unsigned int unique_count = 0;
+    for (unsigned int i = 0; i < *size; i++) {
+        unsigned int is_duplicate = 0;
+        // 检查当前元素是否已经在 unique_indices 中
+        for (unsigned int j = 0; j < unique_count; j++) {
+            if ((*indices)[i] == unique_indices[j]) {
+                is_duplicate = 1;
+                break;
+            }
+        }
+        // 如果没有重复，添加到 unique_indices
+        if (!is_duplicate) {
+            unique_indices[unique_count] = (*indices)[i];
+            unique_count++;
+        }
+    }
+    // 释放原有的 indices
+    free(*indices);
+    // 更新 indices 指针和大小
+    *indices = unique_indices;
+    *size = unique_count;
+}
+
+bool line_vacuate_c(const studio_line_c *line, const int max_points, const double epsilon, studio_line_c *vacuate_line)
+{
+    // 初始化状态变量
+    bool status = false;
+
+    // 如果输入的线段为空，则打印错误信息并返回false
+    if (line->size == 0)
+    {
+        printf("line is empty\n");
+        return false;
+    }
+
+    // 如果输入的线段点数小于等于最大点数，则直接将线段点添加到vacuate_line中，并返回true
+    if (line->size <= max_points)
+    {
+        for (size_t i = 0; i < line->size; ++i)
+        {
+            studio_line_c_add_point(vacuate_line, line->data[i]);
+        }
+        status = true;
+        return status;
+    }
+
+    // 分配内存空间存储线段点的索引
+    unsigned int *indices = (unsigned int *)malloc(line->size * sizeof(unsigned int));
+    unsigned int index_count = 0;
+
+    double t_epsilon = epsilon;
+
+    // 循环执行Douglas-Peucker算法，直到索引点数小于等于最大点数
+    while (1)
+    {
+        index_count = 0;
+        douglas_peucker_c(line, 0, line->size - 1, t_epsilon, indices, &index_count);
+        // 将indices 中的数据去重
+        remove_duplicates(&indices, &index_count);
+
+        if (index_count <= (unsigned int)max_points)
+        {
+            break;
+        }
+        t_epsilon *= 1.1;
+    }
+
+    // 将索引点添加到vacuate_line中
+    for (unsigned int i = 0; i < index_count; ++i)
+    {
+        studio_line_c_add_point(vacuate_line, line->data[indices[i]]);
+    }
+
+    // 释放内存空间
+    free(indices);
+    status = true;
+    return status;
+}

line_vacuate/studio_geo_algo_c.h

@@ -0,0 +1,25 @@
+/**
+ ******************************************************************************
+ * @file           : studio_geo_algo_c.h
+ * @author         : wangyingjie
+ * @brief          : C 语言版 几何算法库： DP算法
+ * @attention      : None
+ * @date           : 2025/5/11
+ ******************************************************************************
+ */
+
+#ifndef STUDIO_GEO_ALGO_C_H
+#define STUDIO_GEO_ALGO_C_H
+
+#include "studio_geo_c.h"
+
+/////////////////// 矢量(线)抽稀算法 ///////////////////
+
+double point_line_dist_square_c(const studio_point_c *p, const studio_point_c *start, const studio_point_c *end);
+
+// Douglas-Peucker算法，用于简化多边形或曲线
+void douglas_peucker_c(const studio_line_c *points, size_t start, size_t end, double epsilon, unsigned int *indices, unsigned int *index_count);
+
+bool line_vacuate_c(const studio_line_c *line, const int max_points, const double epsilon, studio_line_c *vacuate_line);
+
+#endif  // STUDIO_GEO_ALGO_C_H

line_vacuate/studio_geo_c.c

@@ -0,0 +1,167 @@
+/**
+******************************************************************************
+* @file           : studio_geo_c.h
+* @author         : wangyingjie
+* @brief          : 几何矢量类型 C 语言版
+* @attention      : None
+* @date           : 2025/5/10
+******************************************************************************
+*/
+
+#include "studio_geo_c.h"
+
+////////////// 点 //////////////
+
+studio_point_c studio_point_init(double xx, double yy)
+{
+    studio_point_c point = {xx, yy};
+    return point;
+}
+
+bool studio_point_equal(const studio_point_c *p1, const studio_point_c *p2)
+{
+    bool is_equal = fabs(p1->x - p2->x) <= AO_EPSILON && fabs(p1->y - p2->y) <= AO_EPSILON;
+    return is_equal;
+}
+
+////////////// 线 //////////////
+
+// 初始化线段 (默认容量为4)
+studio_line_c studio_line_c_init()
+{
+    studio_line_c line;
+    line.size = 0;
+    line.capacity = 0;
+    line.data = (studio_point_c *)malloc(line.capacity * sizeof(studio_point_c));
+    return line;  // 注意: 返回结构体副本 (需配合 C99 或更高标准)
+}
+
+// 销毁线段 (释放内存)
+void studio_line_c_destroy(studio_line_c *line)
+{
+    if (line->data)
+    {
+        free(line->data);
+        line->data = NULL;
+        line->size = 0;
+        line->capacity = 0;
+    }
+}
+
+// 添加点 (自动扩容)
+void studio_line_c_add_point(studio_line_c *line, studio_point_c point)
+{
+    // 容量不足时扩容 (策略: 增加固定容量 10)
+    if (line->size >= line->capacity)
+    {
+        line->capacity += 10;
+        studio_point_c *new_data = (studio_point_c *)realloc(line->data, line->capacity * sizeof(studio_point_c));
+        if (!new_data)
+        {
+            // 此处可添加错误处理 (例如: 终止程序或返回错误码)
+            return;
+        }
+        line->data = new_data;
+    }
+
+    // 添加新元素
+    line->data[line->size++] = point;
+}
+
+// 获取当前元素数量
+unsigned int studio_line_c_size(const studio_line_c *line)
+{
+    return line->size;
+}
+
+// 获取指定位置的点的引用 (注意索引越界问题)
+const studio_point_c *studio_line_c_get_point(const studio_line_c *line, unsigned int index)
+{
+    if (index < line->size)
+    {
+        return &(line->data[index]);
+    }
+    return NULL;  // 越界返回空指针
+}
+
+////////////// 矩形 //////////////
+
+studio_rect_c studio_rect_init(double l, double t, double r, double b)
+{
+    studio_rect_c rect = {studio_point_init(l, t), studio_point_init(r, b)};
+    return rect;
+}
+
+void studio_rect_correct(studio_rect_c *rect)
+{
+    if (rect->left_top.x > rect->right_bottom.x)
+    {
+        double temp = rect->left_top.x;
+        rect->left_top.x = rect->right_bottom.x;
+        rect->right_bottom.x = temp;
+    }
+    if (rect->left_top.y < rect->right_bottom.y)
+    {
+        double temp = rect->left_top.y;
+        rect->left_top.y = rect->right_bottom.y;
+        rect->right_bottom.y = temp;
+    }
+}
+
+bool studio_rect_intersect(const studio_rect_c *r1, const studio_rect_c *r2)
+{
+    return !(r1->right_bottom.x < r2->left_top.x || r1->left_top.x > r2->right_bottom.x || r1->right_bottom.y > r2->left_top.y || r1->left_top.y < r2->right_bottom.y);
+}
+
+////////////// 圆 //////////////
+
+studio_circle_c studio_circle_init(studio_point_c center, double radius)
+{
+    studio_circle_c circle = {center, radius};
+    return circle;
+}
+
+double studio_circle_area(const studio_circle_c *circle)
+{
+    return AO_M_PI * circle->radius * circle->radius;
+}
+
+////////////// 三角形 //////////////
+
+studio_triangle_c studio_triangle_init(studio_point_c a, studio_point_c b, studio_point_c c)
+{
+    studio_triangle_c triangle = {a, b, c};
+    return triangle;
+}
+
+double studio_triangle_oriented_area(const studio_triangle_c *triangle)
+{
+    return (triangle->a.x * (triangle->b.y - triangle->c.y) + triangle->b.x * (triangle->c.y - triangle->a.y) + triangle->c.x * (triangle->a.y - triangle->b.y)) / 2.0;
+}
+
+double studio_triangle_area(const studio_triangle_c *triangle)
+{
+    return fabs(studio_triangle_oriented_area(triangle));
+}
+
+////////////// 椭圆 //////////////
+
+studio_ellipse_c studio_ellipse_init(studio_point_c center, double rx, double ry)
+{
+    studio_ellipse_c ellipse = {center, rx, ry};
+    return ellipse;
+}
+
+double studio_ellipse_area(const studio_ellipse_c *ellipse)
+{
+    return AO_M_PI * ellipse->rx * ellipse->ry;
+}
+
+double studio_ellipse_circumference(const studio_ellipse_c *ellipse, int Ramanujan)
+{
+    if (Ramanujan == 1)
+    {
+        return AO_M_PI * (3 * (ellipse->rx + ellipse->ry) - sqrt((3 * ellipse->rx + ellipse->ry) * (ellipse->rx + 3 * ellipse->ry)));
+    }
+    return 2 * AO_M_PI * sqrt((ellipse->rx * ellipse->rx + ellipse->ry * ellipse->ry) / 2);
+}

line_vacuate/studio_geo_c.h

@@ -0,0 +1,161 @@
+/**
+ ******************************************************************************
+ * @file           : studio_geo_c.h
+ * @author         : wangyingjie
+ * @brief          : 几何矢量类型 C 语言版
+ *                  包括： 点、多点、线、环、面 、圆、 三角、 椭圆
+ * @attention      : None
+ * @date           : 2025/5/10
+ ******************************************************************************
+ */
+
+#ifndef STUDIO_GEO_C_H
+#define STUDIO_GEO_C_H
+
+#include <stdio.h>
+#include <math.h>
+#include <stdlib.h>
+#include <stdbool.h>
+
+#define AO_EPSILON 1e-6
+#define AO_M_PI 3.14159265358979323846
+
+/****************************************/
+/// 点
+/****************************************/
+typedef struct
+{
+    double x;
+    double y;
+} studio_point_c;
+
+/// 初始化点
+/// \param xx
+/// \param yy
+/// \return
+studio_point_c studio_point_init(double xx, double yy);
+
+/// 判断两个点是否相等
+/// \param p1
+/// \param p2
+/// \return true 表示 相等
+bool studio_point_equal(const studio_point_c *p1, const studio_point_c *p2);
+
+/****************************************/
+/// 线
+/****************************************/
+
+typedef struct
+{
+    studio_point_c *data;   // 点数据指针
+    unsigned int size;      // 当前元素个数
+    unsigned int capacity;  // 当前分配的内存容量
+} studio_line_c;
+
+// 初始化线段 (默认容量为4)
+studio_line_c studio_line_c_init();
+
+// 销毁线段 (释放内存)
+void studio_line_c_destroy(studio_line_c *line);
+
+// 添加点 (自动扩容)
+void studio_line_c_add_point(studio_line_c *line, studio_point_c point);
+
+// 获取当前元素数量
+unsigned int studio_line_c_size(const studio_line_c *line);
+
+// 获取指定位置的点的引用 (注意索引越界问题)
+const studio_point_c *studio_line_c_get_point(const studio_line_c *line, unsigned int index);
+
+/****************************************/
+/// 矩形框
+/****************************************/
+typedef struct
+{
+    studio_point_c left_top;
+    studio_point_c right_bottom;
+} studio_rect_c;
+
+studio_rect_c studio_rect_init(double l, double t, double r, double b);
+
+/// 矩形框校正
+/// \param rect
+void studio_rect_correct(studio_rect_c *rect);
+
+/// 矩形框是否相交
+/// \param r1
+/// \param r2
+/// \return
+bool studio_rect_intersect(const studio_rect_c *r1, const studio_rect_c *r2);
+
+/****************************************/
+/// 圆
+/****************************************/
+typedef struct
+{
+    studio_point_c center;
+    double radius;
+} studio_circle_c;
+
+/// 初始化圆
+/// \param center 圆心
+/// \param radius 半径
+/// \return
+studio_circle_c studio_circle_init(studio_point_c center, double radius);
+
+/// 计算圆的面积
+/// \param circle
+/// \return
+double studio_circle_area(const studio_circle_c *circle);
+
+/****************************************/
+/// 三角形
+/****************************************/
+typedef struct
+{
+    studio_point_c a;
+    studio_point_c b;
+    studio_point_c c;
+} studio_triangle_c;
+
+/// 初始化三角形
+/// \param a
+/// \param b
+/// \param c
+/// \return
+studio_triangle_c studio_triangle_init(studio_point_c a, studio_point_c b, studio_point_c c);
+
+/// 计算三角形的面积
+/// \param triangle
+/// \return
+double studio_triangle_oriented_area(const studio_triangle_c *triangle);
+
+/// 计算三角形的面积
+/// \param triangle
+/// \return
+double studio_triangle_area(const studio_triangle_c *triangle);
+
+/****************************************/
+/// 椭圆
+/****************************************/
+typedef struct
+{
+    studio_point_c center;
+    double rx;
+    double ry;
+} studio_ellipse_c;
+
+studio_ellipse_c studio_ellipse_init(studio_point_c center, double rx, double ry);
+
+/// 计算椭圆的面积
+/// \param ellipse
+/// \return
+double studio_ellipse_area(const studio_ellipse_c *ellipse);
+
+/// 计算椭圆的周长
+/// \param ellipse
+/// \param Ramanujan
+/// \return
+double studio_ellipse_circumference(const studio_ellipse_c *ellipse, int Ramanujan);
+
+#endif  //  STUDIO_GEO_C_H

line_vacuate/task_algo_c.cpp

@@ -0,0 +1,121 @@
+/**
+******************************************************************************
+* @file           : task_alog_c.cpp
+* @author         : wyj
+* @brief          : C语言语法测试
+* @attention      : None
+* @date           : 2025/5/9
+******************************************************************************
+*/
+
+#ifdef __cplusplus
+extern "C"
+{
+#endif
+
+#include "geography/studio_proj_c.h"
+#include "geometry/studio_geo_algo_c.h"
+
+#ifdef __cplusplus
+}
+#endif
+
+#include <stdio.h>
+#include <stdio.h>
+#include "geometry/studio_geo_utils.h"
+
+void studio_line_to_c_line(const studio_line &line, studio_line_c *c_line)
+{
+//    c_line->size = line.size();
+    c_line->data = (studio_point_c *)malloc(line.size() * sizeof(studio_point_c));
+    for (size_t i = 0; i < line.size(); ++i)
+    {
+        studio_point_c p = studio_point_init(line[i].x, line[i].y);
+        studio_line_c_add_point(c_line, p);
+    }
+}
+
+void c_line_to_studio_line(const studio_line_c *c_line, studio_line &line)
+{
+    for (size_t i = 0; i < c_line->size; ++i)
+    {
+        const studio_point_c *p = studio_line_c_get_point(c_line, i);
+        line.push_back(studio_point(p->x, p->y));
+    }
+}
+
+int main()
+{
+    printf("\n\n===================== %s =====================\n\n", __FILE__);
+    silly::geo::utils::init_gdal_env();
+    std::string path;
+#ifdef IS_WINDOWS
+    path = "D:/5_file/2_readfile/geojson/multi_point/fitting_examples.geojson";
+#else
+    path = "/home/wyj/myself/2_data/2_geojson/multi_point/fitting_examples.geojson";
+#endif
+
+    std::cout << "path: " << path << std::endl;
+    std::vector<studio_geo_coll> res_collections;
+    std::vector<studio_geo_coll> collections;
+    silly::geo::utils::read_geo_coll(path, collections);
+    for (auto &coll : collections)
+    {
+        // ------------- 转换为高斯投影 -------------
+        studio_line gauss_line;
+        double central = static_cast<int>(coll.m_line[0].x / 3) * 3;
+        for (auto &point : coll.m_line)
+        {
+            double gx = 0.0;
+            double gy = 0.0;
+            lonlat_to_gauss(central, point.x, point.y, &gx, &gy);
+            gauss_line.push_back(studio_point(gx, gy));
+        }
+        studio_line_c gauss_line_c = studio_line_c_init();
+        studio_line_c vac_gauss_line_c = studio_line_c_init();
+        studio_line_to_c_line(gauss_line, &gauss_line_c);
+
+        // 简化线段，目标点数为28个
+        int max_points = 28;
+        double epsilon = 1.0;
+
+        bool res = line_vacuate_c(&gauss_line_c, max_points, epsilon, &vac_gauss_line_c);
+        if (!res)
+        {
+            std::cout << "Failed to simplify line." << std::endl;
+            return 1;
+        }
+
+        // 高斯投影在转回经纬度
+        studio_line simplified_line;
+        for (int i = 0; i < vac_gauss_line_c.size; i++)
+        {
+            double lon = 0.0;
+            double lat = 0.0;
+            const studio_point_c *p = studio_line_c_get_point(&vac_gauss_line_c, i);
+            gauss_to_lonlat(central, p->x, p->y, &lon, &lat);
+            simplified_line.push_back(studio_point(lon, lat));
+        }
+
+        studio_geo_coll temp;
+        temp.m_type = enum_geometry_type::egtLineString;
+        temp.m_line = simplified_line;
+        res_collections.push_back(temp);
+
+        studio_line_c_destroy(&gauss_line_c);
+        studio_line_c_destroy(&vac_gauss_line_c);
+
+        break;
+    }
+    std::string output_path;
+#ifdef IS_WINDOWS
+    output_path = "D:/5_file/2_readfile/geojson/multi_point/fitting_examples_res_1.geojson";
+#else
+    output_path = "/home/wyj/myself/2_data/2_geojson/multi_point/fitting_examples_res_1_c.geojson";
+#endif
+    silly::geo::utils::write_geo_coll(output_path, res_collections);
+
+    silly::geo::utils::destroy_gdal_env();
+
+    return 0;
+}