Obtaining Least square adjusted single line by intersecting many 3D planes

Question

I am working with many 3D planes and looking for a Least square solution for below case.

IF I am having many number of 3D planes knowing only one point and the normal vector (for eg. O1 and N1), and all those planes intersect each other and make nearly very close 3d lines, then how to compute the least square adjusted one single 3d line to represent all those intersections.

To get a clear idea, I have inserted a figure.

• Known: one point and the normal vector of each plane.
• Find: Least Square fitted single line3d

as i want to do this with c++, i use c++ tag too.

Answer 1

Entirely untested.

What if you took the directions of the lines from the intersections and got the Principle Component

This would get you the direction they're going in. Then create a plane using that direction and an arbitrary point, Project all the points from the plane intersection calculations onto the plane, and find the mean point of these projected points.

Use that mean point and the principle component to define your line.

Something like...

class Plane
{
public:
    Vector3 Point;
    Vector3 Normal;

    Line Intersect (const Plane &other);

    Vector3 Project (const Vector3 &point);
}

class Line
{
public:
    Vector3 Point;
    Vector3 Direction;

    Line (Vector3 point, Vector3 dir);

};

Vector3 PrincipleComponent (const std::vector<Line> &lines)
{
    //You could use the covariance matrix to get this but I will try the interative method on wikipedia.
    Vector3 p(1,2,3); //a random vector?
    static const int c = 10;
    for (int i = 0; i < c; ++i)
    {
        Vector3 t;
        for (auto i = lines.begin(); i != lines.end (); ++i)
        {
            t = t + ((*i).Direction.Dot (p)) * (*i).Direction;
        }
        t.Normalize();
        p = t;
    }
    return p;
}

int main ()
{
    std::vector<Line> LinesFromPlaneIntersections;


    Vector3 direction = PrincipleComponent (LinesFromPlaneIntersections);
    Plane projplane;
    projplane.Normal = direction;
    projplane.Point = LinesFromPlaneIntersections[0].Point;

    Vector3 meanpoint;
    for (auto i = LinesFromPlaneIntersections.begin(); i != LinesFromPlaneIntersections.end (); ++i)
    {
        meanpoint += projplane.Project ((*i).Point);
    }

    meanpoint /= LinesFromPlaneIntersections.size ();

    Line result (meanpoint,direction);
}