背景

前几天写了一篇关于旋转门的数据压缩算法在PostgreSQL中的实现，里面用到了PostGIS里面的ST_Azimuth函数用来计算夹角，其实在PostgreSQL 中，我们还可以使用三角函数，以及三边来求夹角。

文中用到的计算夹角的方法如下

SELECT 180-ST_Azimuth(                              ST_MakePoint(o_x, o_val+i_radius),    -- 门上点                              ST_MakePoint(v_x, v_val)              -- next point                           )/(2*pi())*360 as degAz,                 -- 上夹角                 ST_Azimuth(                              ST_MakePoint(o_x, o_val-i_radius),    -- 门下点                              ST_MakePoint(v_x, v_val)              -- next point                           )/(2*pi())*360 As degAzrev               -- 下夹角      INTO v_angle1, v_angle2;

余弦定理

cosA=(b²+c²-a²)/(2bc)

定点为A、B、C； 对的边分别为a、b、c；

PostgreSQL 支持的三角函数

例子

已知三个点A(3,2)，B(1,2.5)，C(1,1)。 求夹角B, C。

套用余弦公式

cosB=(a²+c²-b²)/(2ac)     cosC=(b²+a²-c²)/(2ba)

首先求三条边长

postgres=# select point_distance(point(3,2), point(1,2.5)) as c , point_distance(point(3,2), point(1,1)) as b , point_distance(point(1,1), point(1,2.5)) as a;        c         |        b         |  a  ------------------+------------------+----- 2.06155281280883 | 2.23606797749979 | 1.5(1 row)

运算如下

cosB=(a²+c²-b²)/(2ac)      =(1.5^2 + 2.06155281280883^2 - 2.23606797749979^2) / (2*1.5*2.06155281280883)  =0.24253562503633260164    cosC=(b²+a²-c²)/(2ba)     =(1.5^2 + 2.23606797749979^2 - 2.06155281280883^2) / (2*2.23606797749979*1.5)  =0.44721359549995825124

求夹角 1 度数

postgres=# select acosd(0.24253562503633260164);      acosd       ------------------ 75.9637565320735(1 row)

求夹角 2 度数

postgres=# select acosd(0.44721359549995825124);      acosd      ----------------- 63.434948822922(1 row)

比对使用PostGIS计算的结果一致

test=>  SELECT 180-ST_Azimuth(                              ST_MakePoint(1,2.5),    -- 门上点                              ST_MakePoint(3,2)              -- next point                           )/(2*pi())*360 as degAz,                 -- 上夹角                 ST_Azimuth(                              ST_MakePoint(1,1),    -- 门下点                              ST_MakePoint(3,2)              -- next point                           )/(2*pi())*360 As degAzrev ;      degaz       |    degazrev     ------------------+----------------- 75.9637565320735 | 63.434948822922(1 row)

源码

三角函数属于浮点运算中的函数

src/backend/utils/adt/float.c

/* *              acosd_q1                - returns the inverse cosine of x in degrees, for x in *                                                the range [0, 1].  The result is an angle in the *                                                first quadrant --- [0, 90] degrees. * *                                                For the 3 special case inputs (0, 0.5 and 1), this *                                                function will return exact values (0, 60 and 90 *                                                degrees respectively). */static doubleacosd_q1(double x){        /*         * Stitch together inverse sine and cosine functions for the ranges [0,         * 0.5] and (0.5, 1].  Each expression below is guaranteed to return         * exactly 60 for x=0.5, so the result is a continuous monotonic function         * over the full range.         */        if (x <= 0.5)        {                volatile float8 asin_x = asin(x);                return 90.0 - (asin_x / asin_0_5) * 30.0;        }        else        {                volatile float8 acos_x = acos(x);                return (acos_x / acos_0_5) * 60.0;        }}/* *              dacosd                  - returns the arccos of arg1 (degrees) */Datumdacosd(PG_FUNCTION_ARGS){        float8          arg1 = PG_GETARG_FLOAT8(0);        float8          result;        /* Per the POSIX spec, return NaN if the input is NaN */        if (isnan(arg1))                PG_RETURN_FLOAT8(get_float8_nan());        INIT_DEGREE_CONSTANTS();        /*         * The principal branch of the inverse cosine function maps values in the         * range [-1, 1] to values in the range [0, 180], so we should reject any         * inputs outside that range and the result will always be finite.         */        if (arg1 < -1.0 || arg1 > 1.0)                ereport(ERROR,                                (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),                                 errmsg("input is out of range")));        if (arg1 >= 0.0)                result = acosd_q1(arg1);        else                result = 90.0 + asind_q1(-arg1);        CHECKFLOATVAL(result, false, true);        PG_RETURN_FLOAT8(result);}

man asin NAME       asin, asinf, asinl - arc sine functionSYNOPSIS       #include 
    
            double asin(double x);       float asinf(float x);       long double asinl(long double x);///CONFORMING TO       C99, POSIX.1-2001.  The variant returning double also conforms to SVr4, 4.3BSD, C89.SEE ALSO       acos(3), atan(3), atan2(3), casin(3), cos(3), sin(3), tan(3)
    

祝大家玩得开心，欢迎随时来 阿里云促膝长谈业务需求 ，恭候光临。

阿里云的小伙伴们加油，努力 做好内核与服务，打造最贴地气的云数据库 。