空间直线同球体交点求解

 

 一、问题求解

空间直线同球体交点问题求解，为求解二元二次方程问题，下面是求解过程

 

二、代码实现

1、代码实现

void cal_q ( double par_x, double par_y, double par_z, double par_r, double bpx, double bpy, double bpz, double spx, double spy, double spz ) {
    double x_1, x_2, y_1, y_2, z_1, z_2;

    if ( bpx != spx ) {
        double k_xy = ( spy - bpy ) / ( spx - bpx );
        double b_xy = bpy - k_xy * bpx;

        double k_zx = ( spz - bpz ) / ( spx - bpx );
        double b_zx = bpz - k_zx * bpx;

        double A = 1 + k_xy * k_xy + k_zx * k_zx;
        double B = 2 * k_xy * ( b_xy - par_y ) + 2 * k_zx * ( b_zx - par_z ) - 2 * par_x;
        double C = par_x * par_x + ( b_xy - par_y ) * ( b_xy - par_y ) + ( b_zx - par_z ) * ( b_zx - par_z ) - par_r * par_r;

        x_1 = ( -B + sqrt ( B * B - 4 * A * C ) ) / ( 2 * A );
        y_1 = k_xy * x_1 + b_xy;
        z_1 = k_zx * x_1 + b_zx;

        x_2 = ( -B - sqrt ( B * B - 4 * A * C ) ) / ( 2 * A );
        y_2 = k_xy * x_2 + b_xy;
        z_2 = k_zx * x_2 + b_zx;

        printf ( "x_1: %10.5f   y_1: %10.5f   z_1: %10.5f \n", x_1, y_1, z_1 );
        printf ( "x_2: %10.5f   y_2: %10.5f   z_2: %10.5f \n", x_2, y_2, z_2 );

        printf ( "bpx != spx \n" );
    } else {
        if ( bpy != spy ) {
            double k_zy = ( spz - bpz ) / ( spy - bpy );
            double b_zy = spz - k_zy * spy;

            double A = 1 + k_zy * k_zy;
            double B = 2 * k_zy * ( b_zy - par_z ) - 2 * par_y;
            double C = par_y * par_y + ( b_zy - par_z ) * ( b_zy - par_z ) + ( bpx - par_x ) * ( bpx - par_x ) - par_r * par_r;

            x_1 = bpx;
            y_1 = ( -B + sqrt ( B * B - 4 * A * C ) ) / ( 2 * A );
            z_1 = k_zy * y_1 + b_zy;

            x_2 = bpx;
            y_2 = ( -B - sqrt ( B * B - 4 * A * C ) ) / ( 2 * A );
            z_2 = k_zy * y_2 + b_zy;

            printf ( "x_1: %10.5f   y_1: %10.5f   z_1: %10.5f \n", x_1, y_1, z_1 );
            printf ( "x_2: %10.5f   y_2: %10.5f   z_2: %10.5f \n", x_2, y_2, z_2 );

            printf ( "bpx == spx && bpy != spy \n" );
        } else {
            if ( bpz != spz ) {
                x_1 = bpx;
                y_1 = bpy;
                z_1 = par_z + sqrt ( par_r * par_r - ( bpx - par_x ) * ( bpx - par_x ) - ( bpy - par_y ) * ( bpy - par_y ) );

                x_2 = bpx;
                y_2 = bpy;
                z_2 = par_z - sqrt ( par_r * par_r - ( bpx - par_x ) * ( bpx - par_x ) - ( bpy - par_y ) * ( bpy - par_y ) );

                printf ( "x_1: %10.5f   y_1: %10.5f   z_1: %10.5f \n", x_1, y_1, z_1 );
                printf ( "x_2: %10.5f   y_2: %10.5f   z_2: %10.5f \n", x_2, y_2, z_2 );

                printf ( "bpx == spx && bpy == spy \n" );
            } else {
                printf ( "Warning: two points overlap !!!! \n" );64             }
        }
    }
}

2、测试

#include <cstdio>
#include <cmath>

void cal_q ( double, double, double, double, double, double, double, double, double, double );

int main() {
    double par_x = 20;
    double par_y = 20;
    double par_z = 20;
    double par_r = 15;

    double bpx = 20;
    double bpy = 30;
    double bpz = 40;

    double spx = 20;
    double spy = 20;
    double spz = 20;

    printf ( "bpx: %10.5f   bpy: %10.5f   bpz: %10.5f \n", bpx, bpy, bpz );
    printf ( "spx: %10.5f   spy: %10.5f   spz: %10.5f \n", spx, spy, spz );

    cal_q ( par_x, par_y, par_z, par_r, bpx, bpy, bpz, spx, spy, spz );

    return 0;
}

3、图形显示测试结果