arc_From_3pt.c

float *lerp(float *a, float *b, float t){
    if (t == 0.0) { return a; }
    else if (t == 1.0){  return b; }
    float s = 1.0 - t;
    float x = s * a[0] + t * b[0];
    float y = s * a[1] + t * b[1];
    float z = s * a[2] + t * b[2];
    float *outvalue = {x, y, z};
    return outvalue;
}

float *normal(float *p1, float *p2, float *p3){
    float x = ((p2[1]-p1[1])*(p3[2]-p1[2]))-((p2[2]-p1[2])*(p3[1]-p1[1]));
    float y = ((p2[2]-p1[2])*(p3[0]-p1[0]))-((p2[0]-p1[0])*(p3[2]-p1[2]));
    float z = ((p2[0]-p1[0])*(p3[1]-p1[1]))-((p2[1]-p1[1])*(p3[0]-p1[0]));
    float *outvalue = {x, y, z};
    return outvalue;    
}

float *sub2v3(float *p1, float *p2){
    float *outvalue = {p1[0]-p2[0], p1[1]-p2[1], p1[2]-p2[2]};
    return outvalue;
}

float *add2v3(float *p1, float *p2){
    float *outvalue = {p1[0]+p2[0], p1[1]+p2[1], p1[2]+p2[2]};
    return outvalue;
}

float *matrix_multiply(float *v, float M[4][4]){
    float x = M[0][0]*v[0] + M[0][1]*v[1] + M[0][2]*v[2] + M[0][3];
    float y = M[1][0]*v[0] + M[1][1]*v[1] + M[1][2]*v[2] + M[1][3]; 
    float z = M[2][0]*v[0] + M[2][1]*v[1] + M[2][2]*v[2] + M[2][3];
    float *outvalue = {x, y, z};
    return outvalue;
}

void rotationMatrix(float angle, float *axis, float mat_rot[][4]){
    // ...
    // unimplemented
    // ...
    float mat[4][4] = {
        {1.0, 0.0, 0.0, 0.0},
        {0.0, 1.0, 0.0, 0.0},
        {0.0, 0.0, 1.0, 0.0},
        {0.0, 0.0, 0.0, 1.0}
    };
}

float angle2v3(float *a, float *b, float *fallback){
    
    return fallback;
}

void intersect_line_line(float *v1_, float *v2_, float *v3_, float *v4_, float result[2][3]){
    result = {NULL, NULL};

}

void generate3PT(float[] *pts, int num_verts, float *verts){
    /*
    Arc from [start - through - end]
    - call this function only if you have 3 pts,
    - do your error checking before passing to it.
    */
    num_verts -= 1;

    // construction
    float *v1 = pts[0];
    float *v2 = pts[1];
    float *v3 = pts[1];
    float *v4 = pts[2];
    float *edge1_mid = lerp(v1, v2, 0.5);
    float *edge2_mid = lerp(v3, v4, 0.5);
    float *axis = normal(v1, v2, v4);

    float mat_rot[4][4];
    rotationMatrix(math.radians(90.0), axis, mat_rot);

    // triangle edges
    float *v1_ = add2v3(matrix_multiply(sub2v3(v1, edge1_mid), mat_rot), edge1_mid);
    float *v2_ = add2v3(matrix_multiply(sub2v3(v2, edge1_mid), mat_rot), edge1_mid);
    float *v3_ = add2v3(matrix_multiply(sub2v3(v3, edge2_mid), mat_rot), edge2_mid);
    float *v4_ = add2v3(matrix_multiply(sub2v3(v4, edge2_mid), mat_rot), edge2_mid);

    float result[2][3] = {NULL, NULL};
    intersect_line_line(v1_, v2_, v3_, v4_, result);
    
    if (result[0] != NULL) {
        // do arc
        float *p1 = result[0];

        // find arc angle.
        float a = angle2v3(sub2v3(v1, p1), sub2v3(v4, p1), 0);
        float s = (2 * math.pi) - a;

        float interior_angle = angle2v3(sub2v3(v1, v2), sub2v3(v4, v3), 0);
        if (interior_angle > 0.5 * math.pi) {
            s = math.pi + 2 * (0.5 * math.pi - interior_angle);
        }
        
        float *vec;
        float mat_rot2[4][4];
        int i;
        for (i = 0; i < (num_verts + 1); i = i + 1) {
            rotationMatrix(((s / num_verts) * i), axis, mat_rot2);
            vec = add2v3(matrix_multiply(sub2v3(v4, p1), mat_rot), p1);
            verts[i] = vec;
        }

    } else {
    
        // do straight line
        int i;
        float *vec;
        float step_size = 1.0 / float(num_verts);
        for (i = 0; i < (num_verts + 1); i = i + 1) {
            vec = lerp(v1_, v4_, i * step_size);
            verts[i] = vec;
        }

int main () {

    int num_verts = 20;
    float *verts[num_verts-1][3];
    generate3PT(pts, num_verts, verts);

    printf(verts);

    return 0;
};