include/tinyhtm/geometry.h

Go to the documentation of this file.

#ifndef HTM_GEOMETRY_H
#define HTM_GEOMETRY_H

#include "common.h"

#ifdef __cplusplus
extern "C" {
#endif


/* ================================================================ */
/* ================================================================ */

struct htm_v3 {
    double x; 
    double y; 
    double z; 
};

struct htm_sc {
    double lon; 
    double lat; 
};

HTM_INLINE enum htm_errcode htm_v3_init(struct htm_v3 *out,
                                        double x,
                                        double y,
                                        double z)
{
    if (HTM_ISSPECIAL(x) || HTM_ISSPECIAL(y) || HTM_ISSPECIAL(z)) {
        return HTM_ENANINF;
    } else if (out == NULL) {
        return HTM_ENULLPTR;
    }
    out->x = x;
    out->y = y;
    out->z = z;
    return HTM_OK;
}

HTM_INLINE enum htm_errcode htm_sc_init(struct htm_sc *out,
                                        double lon_deg,
                                        double lat_deg)
{
    if (HTM_ISSPECIAL(lon_deg) || HTM_ISSPECIAL(lat_deg)) {
        return HTM_ENANINF;
    } else if (lat_deg < -90.0 || lat_deg > 90.0) {
        return HTM_ELAT;
    } else if (out == NULL) {
        return HTM_ENULLPTR;
    }
    out->lon = lon_deg;
    out->lat = lat_deg;
    return HTM_OK;
}

HTM_INLINE void htm_v3_add(struct htm_v3 *out,
                           const struct htm_v3 *v1,
                           const struct htm_v3 *v2)
{
    out->x = v1->x + v2->x;
    out->y = v1->y + v2->y;
    out->z = v1->z + v2->z;
}

HTM_INLINE void htm_v3_sub(struct htm_v3 *out,
                           const struct htm_v3 *v1,
                           const struct htm_v3 *v2)
{
    out->x = v1->x - v2->x;
    out->y = v1->y - v2->y;
    out->z = v1->z - v2->z;
}

HTM_INLINE void htm_v3_neg(struct htm_v3 *out, const struct htm_v3 *v)
{
    out->x = - v->x;
    out->y = - v->y;
    out->z = - v->z;
}

HTM_INLINE void htm_v3_mul(struct htm_v3 *out,
                           const struct htm_v3 *v,
                           double s)
{
    out->x = v->x * s;
    out->y = v->y * s;
    out->z = v->z * s;
}

HTM_INLINE void htm_v3_div(struct htm_v3 *out,
                           const struct htm_v3 *v,
                           double s)
{
    out->x = v->x / s;
    out->y = v->y / s;
    out->z = v->z / s;
}

HTM_INLINE double htm_v3_dot(const struct htm_v3 *v1,
                             const struct htm_v3 *v2)
{
    return v1->x * v2->x + v1->y * v2->y + v1->z * v2->z;
}

HTM_INLINE void htm_v3_cwise_mul(struct htm_v3 *out,
                                 const struct htm_v3 *v1,
                                 const struct htm_v3 *v2)
{
    out->x = v1->x * v2->x;
    out->y = v1->y * v2->y;
    out->z = v1->z * v2->z;
}

HTM_INLINE double htm_v3_norm2(const struct htm_v3 *v)
{
    return v->x * v->x + v->y * v->y + v->z * v->z;
}

HTM_INLINE double htm_v3_norm(const struct htm_v3 *v)
{
    return sqrt(htm_v3_norm2(v));
}

HTM_INLINE void htm_v3_normalize(struct htm_v3 *out, const struct htm_v3 *v)
{
    double norm = htm_v3_norm(v);
    out->x = v->x / norm;
    out->y = v->y / norm;
    out->z = v->z / norm;
}

HTM_INLINE void htm_v3_rcross(struct htm_v3 *out,
                              const struct htm_v3 *v1,
                              const struct htm_v3 *v2)
{
    double x1 = v2->x + v1->x;
    double x2 = v2->x - v1->x;
    double y1 = v2->y + v1->y;
    double y2 = v2->y - v1->y;
    double z1 = v2->z + v1->z;
    double z2 = v2->z - v1->z;
    out->x = y1 * z2 - z1 * y2;
    out->y = z1 * x2 - x1 * z2;
    out->z = x1 * y2 - y1 * x2;
}

HTM_INLINE void htm_v3_cross(struct htm_v3 *out,
                             const struct htm_v3 *v1,
                             const struct htm_v3 *v2)
{
    double x = v1->y * v2->z - v1->z * v2->y;
    double y = v1->z * v2->x - v1->x * v2->z;
    double z = v1->x * v2->y - v1->y * v2->x;
    out->x = x;
    out->y = y;
    out->z = z;
}

enum htm_errcode htm_v3_ne(struct htm_v3 *north,
                           struct htm_v3 *east,
                           const struct htm_v3 *v);

enum htm_errcode htm_v3_tanrot(double *angle,
                               const struct htm_v3 *v1,
                               const struct htm_v3 *v2,
                               double r);

enum htm_errcode htm_v3_rot(struct htm_v3 *out,
                            const struct htm_v3 *v,
                            const struct htm_v3 *k,
                            double angle_deg);

enum htm_errcode htm_v3_centroid(struct htm_v3 *cen,
                                 const struct htm_v3 *points,
                                 size_t n);

enum htm_errcode htm_sc_tov3(struct htm_v3 *out, const struct htm_sc *p);

enum htm_errcode htm_v3_tosc(struct htm_sc *out, const struct htm_v3 *v);


/* ================================================================ */
/* ================================================================ */

double htm_sc_dist2(const struct htm_sc *p1, const struct htm_sc *p2);

double htm_sc_angsep(const struct htm_sc *p1, const struct htm_sc *p2);

HTM_INLINE double htm_v3_dist2(const struct htm_v3 *v1,
                               const struct htm_v3 *v2)
{
    struct htm_v3 delta;
    htm_v3_sub(&delta, v1, v2);
    return htm_v3_norm2(&delta);
}

double htm_v3_angsepu(const struct htm_v3 *v1,
                      const struct htm_v3 *v2);

double htm_v3_angsep(const struct htm_v3 *v1, const struct htm_v3 *v2);

double htm_v3_edgedist2(const struct htm_v3 *v,
                        const struct htm_v3 *v1,
                        const struct htm_v3 *v2,
                        const struct htm_v3 *e);


/* ================================================================ */
/* ================================================================ */

struct htm_s2ellipse {
    struct htm_v3 cen; 
    double xx;  
    double yy;  
    double zz;  
    double xy;  
    double xz;  
    double yz;  
    double a;   
};

enum htm_errcode htm_s2ellipse_init(struct htm_s2ellipse *e,
                                    const struct htm_v3 *f1,
                                    const struct htm_v3 *f2,
                                    double a);

enum htm_errcode htm_s2ellipse_init2(struct htm_s2ellipse *e,
                                     const struct htm_v3 *cen,
                                     double a,
                                     double b,
                                     double angle);

HTM_INLINE int htm_s2ellipse_cv3(const struct htm_s2ellipse *e,
                                 const struct htm_v3 *v)
{
    double qf = e->xx * v->x * v->x +
                e->yy * v->y * v->y +
                e->zz * v->z * v->z +
          2.0 * e->xy * v->x * v->y +
          2.0 * e->xz * v->x * v->z +
          2.0 * e->yz * v->y * v->z;
    double dp = htm_v3_dot(&e->cen, v);
    if (e->a <= 90.0) {
        return dp >= 0.0 && qf <= 0.0;
    }
    return dp >= 0.0 || qf >= 0.0;
}


/* ================================================================ */
/* ================================================================ */

struct htm_s2cpoly {
    size_t n;               
    struct htm_v3 vsum;     
    struct htm_v3 ve[];     
};

HTM_INLINE struct htm_v3 * htm_s2cpoly_verts(struct htm_s2cpoly *poly)
{
    return poly->ve;
}

HTM_INLINE struct htm_v3 * htm_s2cpoly_edges(struct htm_s2cpoly *poly)
{
    return poly->ve + poly->n;
}

struct htm_s2cpoly * htm_s2cpoly_init(const struct htm_v3 *verts,
                                      size_t n,
                                      enum htm_errcode *err);

struct htm_s2cpoly * htm_s2cpoly_box(const struct htm_v3 *cen,
                                     double width,
                                     double height,
                                     double angle,
                                     enum htm_errcode *err);

struct htm_s2cpoly * htm_s2cpoly_ngon(const struct htm_v3 *cen,
                                      double r,
                                      size_t n,
                                      enum htm_errcode *err);

struct htm_s2cpoly * htm_s2cpoly_line(const struct htm_v3 *v1,
                                      const struct htm_v3 *v2,
                                      double r,
                                      enum htm_errcode *err);

int htm_s2cpoly_cv3(const struct htm_s2cpoly *poly, const struct htm_v3 *v);

double htm_s2cpoly_area(const struct htm_s2cpoly *poly);

struct htm_s2cpoly * htm_s2cpoly_clone(const struct htm_s2cpoly *poly);

enum htm_errcode htm_s2cpoly_pad(struct htm_s2cpoly *poly, double r);


/* ================================================================ */
/* ================================================================ */

int htm_v3_hemispherical(const struct htm_v3 *points,
                         size_t n,
                         enum htm_errcode *err);

int htm_v3_convex(const struct htm_v3 *points,
                  size_t n,
                  enum htm_errcode *err);

struct htm_s2cpoly * htm_s2cpoly_hull(const struct htm_v3 *points,
                                      size_t n,
                                      enum htm_errcode *err);

#ifdef __cplusplus
}
#endif

#endif /* HTM_GEOMETRY_H */