Skip to content
RebeccaHisey edited this page Oct 16, 2016 · 5 revisions

CISC 330 assignment 1

How to run:

  • load the module corresponding to each set of questions
  • all points/ vectors should be numpy arrays
    • Code for questions 1-4 located in Intersections.py

      • enter the following line into the command line: logic = slicer.modules.IntersectionsWidget.logic
      • Question 1. intersection = logic.LineAndPlane(l0,l1,p0, n)
      • Question 2. (intersection, error) = logic.TwoLines(l0_1,l1_1,l0_2,l0_2)
        • where l0_1, l1_1 are points on line 1 and l0_2 and l1_2 are on line 2
      • Question 3. (intersection, error) = NLines(numLines, listOfLines)
        • where numLines is the number of lines included in the list of lines
        • listOfLines is a list of lines, each line is given by 2 points
      • Question 4. intersections = LineAndEllipsoid(l0, l1, a, b, c)
        • l0, l1 are points on the line
        • a,b,c are the axes of the ellipsoid
    • Code for questions 5-6 located in SimAndRecon.py

      • enter the following line into the command line: logic = slicer.modules.SimAndReconWidget.logic
      • Question 5: (cCenter, cRadius, cError) = SphereRecon(ctr, radius, numPoints, maxOff)
        • ctr and radius are the original center and radius of the sphere
        • numPoints is the number of points to be generated on the sphere
        • maxOff is the maximum offset
        • to only simulate a sphere (not perform the reconstruction) dataPoints = SphereSim(ctr, radius, numPoints, maxOff)
      • Question 6: (holdingPoint, direction, error) = LineRecon(A, B, numPoints, maxOff)
        • A, B are points describing the original line
        • numPoints is the number of points to be generated in the simulation
        • maxOff is the maximum offset
        • to only simulate a sphere (not perform the reconstruction) dataPoints = LineSim(A, B, numPoints, maxOff)
    • Code for questions 7-8 located in Transformations.py

      • enter the following line into the command line: logic = slicer.modules.TransformationsWidget.logic
      • Question 7: (e1,e2,e3,center) = OrthoNormalCoordSystem(A, B, C)
        • A, B, C are points
      • Question 8: transformationMatrix = RigidBodyTransformation(A1, B1, C1, A2, B2, C2)
        • A1,B1,C1 are the original points
        • A2,B2,C2 correspond to A1,B1,C2 after they have been transformed by the transformation matrix

Test cases:

  • Each module contains test cases for each question. These can be viewed by hitting reload and test once the module has been loaded

Methods:

  1. Find the intersection between a line and a plane
  • Given 2 points on a line, start by finding the direction vector and normalizing it by dividing it by it's length
  • The computed intersection is given by the following formula:
    • ci = d*l + l0
    • where ci is the computed intersection
    • l is the unit direction vector of the line
    • l0 is a point on the line
    • d is a scalar
  • d = ((p0 - l) . n )/ (l . n)
    • p0 is a point on the plane
    • l is the unit direction vector of the line
    • n is the normal vector of the plane
  1. Find the intersection of two lines
  • Given 2 points on two different lines, start by finding the direction vector of each line just as in question 1
  • Find a line that is perpendicular to both lines by computing the cross product
  • L1 = P1 + t1*v1
  • L2 = P2 + t2*v2
  • where: L1, P1 are points on the first line
  • L2, P2 are points on the second line
  • t1, t2 are scalars
  • v1, v2 are the unit direction vectors for each line
  • the lines intersect when L1 = L2
  • create the same equation for the vector we got from the cross product of line 1 and line 2
  • set P3 = L2 and L3 = L1 these will be the points that are closest to the intersection
  • perform these substitutions in the equation for the 3rd line
  • solve the system of linear equations to get the t values
  • compute the midpoint between the two lines M = (L1 + L2)/2
  • the error metric I used was the distance from the midpoint to L1
    • since M is halfway between L1 and L2 the distance from M to L2 should be the same
  1. Find the intersection of N lines
  • Find the unit direction vector for each pair of points
  • perform a pairwise comparison to find the intersection point of each pair of lines using the function from question 2
  • compute the average of these points to find the symbolic intersection
  • the error metric is the standard deviation
  • To deal with outliers that could affect the symbolic intersection, check to see if each intersection point between each pair of lines is more than 2 standard deviations away from the symbolic intersection
  • if it is we can eliminate that point and recompute the intersection
  1. Find the intersection of a line with an ellipsoid centred at the origin
  • find the unit direction vector of the line
  • the equation of the ellipsoid is (x^2 / a^2) + (y^2 / b^2) + (z^2 / c^2) = 1
    • need to solve for x, y and z
  • equation for the line is L = P + t*v
  • by expanding the equation we get:
  • ((bcvx)^2 + (acvy)^2 + (abvz)^2)(t^2) + ((2 * (bc)^2 * Pxvx)+(2 * (ac)^2 * Pyvy)+(2 * (ab)^2 Pzvz))t + ((bcPx)^2 + (bcPx)^2 + (bcPx)^2 - (ab*c)^2)
  • we can plug this into the quadratic formula to solve for t
  • if the quadratic formula yields complex roots there is no intersection point
  • if the quadratic formula yields the same value twice there is 1 intersection point
  • if the quadratic formula yields two different values there are two intersection points
  • plug the value(s) for t into the equation of the line to solve for the intersection points
  1. Sphere reconstruction and simulation
  • Method for picking random points on a sphere found at : https://www.jasondavies.com/maps/random-points/
    • generate an evenly spaced list (p) of values between 0.5 and 1
      • should have as many elements as we have points
      • range of 0.5 to 1 means that we will only plot points on the "northern hemisphere" of the sphere
    • randomly the maximum offset for each point between -1* maximum offset and the maximum offset
    • randomly select an angle to rotate around the z axis (phi)
    • for each element P in the list p find angle rotated down from the z axis (psi) = arccos (2*P -1)
    • x coordinate = cos(phi)*sin(psi) + (x coordinate of the center point of the sphere) + offset
    • y coordinate = sin(phi)*sin(psi) + (y coordinate of the center point of the sphere) + offset
    • z coordinate = cos(psi) + (z coordinate of the center point of the sphere) + offset
  • Method for finding the least squares solution for the reconstruction of the sphere found at: http://research.cs.queensu.ca/home/cisc271/pdf/class24.pdf (Dr Ellis, CISC 271 class notes)
    • Create a matrix M whose rows are [-2*xjT 1] where xjT is the row vector of each point in the dataset
    • create a matrix b whose rows are [-1*xj.xj] where xj is a point in the dataset (1 row for each point)
    • solve the linear equations M*[ctr, sigma] = b
      • ctr is the reconstructed center point of the sphere, and sigma is ctr.ctr - radius^2
    • use sigma and the center point to solve for the radius
  • error is the standard deviation between the initial center point and radius and the computed center point and radius
    • The error value increased as the maximum offset increased
  1. Line simulation and reconstruction
  • find the direction vector of the line by subtracting point A - point B and dividing by the length
  • randomly select the offset for each coordinate of each point
  • randomly select a scalar t to change the length of the direction vector by
  • offset point = point A - t*(direction vector) + offset
  • Reconstruction:
    • Used method described in class by Dr. Fichtinger
    • Find the average of all the data points (this is our holding point)
    • compute the vector from the holding point to each of the points in the dataset
      • multiply the vector by -1 if it is pointing in the opposite direction
    • find the average of these vectors, and make it unit length
      • this is the direction vector of the reconstructed line
  • error is the standard deviation between the original vector and the computed vector
    • error stayed mostly consistent as the max offset increased up to 10
  1. Generate an orthonormal coordinate system located at the center of a triangle of points
  • first find the center of gravity by finding the average of the 3 points
  • create 3 vectors from the center to each of the three points
  • the coordinates of these vectors go to form the columns of a matrix
  • perform a QR decomposition of this matrix
  • the columns of Q are vectors that form an orthonormal basis of the vector space
  • tests check that the vectors are orthogonal to each other and all have unit length
  1. Rigid body transformation
  • Method for computing the rotation matrix and translation matrix found at: http://www.comp.nus.edu.sg/~hcheng/academic/HYPprojects/eyeBalls/ch8_image_reg.pdf (page 395)
  • find the average of the 3 points in set 1 (A1,B1, C1) and set 2 (A2, B2 and C2)
  • subtract this average from each of the 3 points in the set
    • the two sets of points are now related by the rotation matrix only
  • use single value decomposition to get the rotation matrix
  • the translation vector is calculated by taking A2 - (rotation matrix)*A1
  • homogeneous translation matrix is 1,0,0,t[0, [0,1,0,t[1]], [0,0,1,t[2]], [0,0,0, 1 ]]
    • where t is the translation vector
  • homogeneous rotation matrix is [[r00, r01, r02, 0], [r10, r11, r12, 0], [r20, r21, r22, 0], [0 , 0 , 0 , 1]]
    • where rij is the value of the ith row and the jth column of the rotation matrix
  • the overall transformation matrix is the homogeneous translation matrix * the homogeneous rotation matrix
Clone this wiki locally