Cells In itkMesh

I have completed inserting the face data into an itkMesh that includes the point data for the ring in the 3D-registration process. By changing the parameters in the registration to use an itkMesh I was able to get the registration successfully running. Then I began to work on the intensity process metric and ran into some problems. Declaring the iterators and different variables in the class using the templated names of the data types has proved difficult. With the inclusion of the itkMesh the iterators for moving through the cells were elusive in the ITK class structure so I had to declare them explicitly. Then I found that there is no simple way to access the cell data in the itkMesh. Some of the reason for this is that ITK allows for a heterogenous sampling of cells to be included in an itkMesh. The methods in the ITK manual include creating visitor classes that seem way too complicated for the task of simply accessing information stored in a cell. I am currently looking for functions that could easily get the point data in the cells so that I will be able to calculate the normals.

3D-Registration

Using a three-dimensional point set requires the registration to be changed so that it contains a three-dimensional transform, optimizer, and image. I first tried to use the itkRigid3DTransform but found that it was also a base class for more specific transforms. itkVersorRigid3DTransform was one of the transforms derived from itkRigid3DTransform and allowed rotation along with translation. The transform had its own tailored optimizer called itkVersorRigid3DTransformOptimizer that would work over the six parameter in the transform. After inserting the optimizer and transform into the registration process I needed to change all the instances of point sets in the registration and visualization into three-dimension ones. The initial parameters in the transform and the scales in the optimizer needed to be augmented to take into account there now being six parameters.

Once the registration was up and running I started to play around with different coordinates to see whether it was working properly. I started by enabling just translation in the x and y-direction by scaling out all the other parameters. Thus the registration was basically two-dimensional over a slice of the three-dimensional image. But the metric value remained extremely small. Throughout the extant of the process there were not any drastic or significant changes. This could be the result of the image not being read correctly into ITK or the point set being in the incorrect position relative to the image.

Displaying The 2D-Slice

When I first tried to display the slice that came from the 3D-image, it resulted in VTK displaying a completely deformed image that had the bulk of the image colored blue with a red border. I had originally thought that this had to do with the itkImageToVTKImageFilter distorting the data received from ITK. Therefore, I wrote out the image as a .png file and loaded that into VTK, but the problem persisted. This led to me to think that the problem stemmed from the format of the data that I was passing to the VTK interface. The 3D-image had been loaded with a pixel type of unsigned short. I had used unsigned char when VTK had previously worked to display .png files. So I went about trying to change the pixel format of the slice that I had removed. The first filter I tried was the itkImageToImageFilter that would not instantiate itself as a new object due to it being a base class to be used only by other classes derived from it. Then I tried to use the itkRescaleIntensityImageFilter filter. This was able to convert from unsigned short pixels to unsigned char and output an image in that format. I connected the output from the filter to the itkImageToVTKImageFilter and passed that into my VTK visualization set up to successfully display the slice.

The following is an image of the distorted slice on the left and the corrected version on the right:

itkPointSet And 3D-Images

The problem with loading the itkPointSet into the Registration process had to do with the point data part of the object being empty. When a point set was outputted from the narrow band filter, it filled the point data for each point with the distance it was from the outside of the band. I am not sure exactly how the registration process uses the data. Although in the case of my ellipse, the registration did not work as well when I inserted large values into the data area of itkPointSet that I was creating. After trying various values, I found the optimal value for the point data to be zero.

To make the ring I need to create a point set with three dimensions. This will require the transform to be changed to something that can handle three dimensional point sets. The image used in the process will also need to be in three dimensions. So I began the process of loading in a three dimensional image. I focused my attention on the visualization part using VTK where a two-dimensional image will be displayed. Getting the slice of the three dimensional image that I wanted was accomplished by using the itkExtractImageFilter. I was able to write the slice out to a .png file, but if I passed it to vtk's visualization machinery using itkImageToVTKImageFilter it produced an image that was completely distorted. To fix this I am trying to pass the two-dimensional slice I have into the itk::ImageToImageFilter to change its pixel representation.

Creating An Ellipse

To create a set of points to represent an ellipse, I first used the equation (x²/79²) + (y²/95²) = 1, solving for y and inserting x-values from -79 to 79, with increments of 0.5. This created around 600 points. The problem was that the sides of the ellipse contained very few points, due to the relatively high changes in the y-direction along the left and right sides of the ellipse. Using the equation in polar coordinates resolved that problem. This was accomplished by incrementing from zero to 2π and plugging that into r² = (79²95²) / (79²sin²(Θ) + 95²cos²(Θ)) as the angle. After converting the resulting polar coordinates into rectangular ones, I was able to create an ellipse that had a uniform spread of points along the border. The following is an image of the ellipse using the rectangular coordinate representation of the equation on the left, and the one using the polar coordinate representation on the right:

Once I had all the points stored into an itkPointSet, I tried to incorporate that into the registration process, but ended up with a segmentation fault error. Taking the itkPointSet that was the output from the spatial object filter on the old ellipse, and replacing its point data with my newly created points fixed the problem. There is some parameter on the itkPointSet that I need to set to get it to work properly. I used the function "Print( std::cout )", that is a part of each itkPointSet, to print out and compare the details of the two itkPointSet objects. The only differences were in the fields "Requested Number Of Regions" with the old point set having its value as two, and mine have a value of zero, and in the field "Reference Count" with a two for the old point set and one for mine.

itkPointSet And itkMesh

I had created an itkPointSetNormal class that was derived from itkPointSet, but now it seems that the functionality I was trying to capture is already in the itkMesh class. My task was to create a point set that contained the normals along with the points. Although calculating the normals would be difficult if the object started with was more complicated than a circle. One solution is to create a ring structure made up of triangular faces that would then be used to calculate the normals of the points. The class itkMesh allows for all of the abilities of itkPointSet with the addition of cells that contain information about the connectivity of two, three, or four points. Each cell can also have a piece of data attached to it. This allows there to be cells with three points representing a face along with a normal attached as the data.

In the current registration, I am using an ellipse as the shape. So to create the ring structure, I will use the equation (x²/79²) + (y²/95²) = 1 which will give me a row of points that represents the boundary of the ellipse, and where all the coordinates in the z-direction are zero. Supplementing that will be a duplicate ellipse except that in this one the z-coordinates will be negative one. This permits the two ellipses to connect together through triangular faces on their perimeter. The calculations of the normals will be handled in the metric class, where the direction will be the average of the normals of the faces that are adjacent to a point.

Calculating Normals

The next step in furthering the registration process is to add intensity processes. Sean Mann, another undergrad in the lab, has already written the retrieval of the intensities from a mesh and image in a class. To incorporate that into the registration I need to pass the point set in as a vtkMesh and calculate the normals from my point set.

I finished writing some functions that calculate the normals from the itkPointSet that I get from using the itkBinaryMaskToNarrowBandPointSetFilter on an itkSpatialObject, which in this case was an ellipse. The functions go through each point and calculate the normal based on the shape being roughly a circle and taking a ray from the center going through a perimeter point. The resulting direction is stored in a float array as a point that is then normalized to a unit vector. This approach seems to work well if the shape is a circle. In the case of other shapes, the normals cannot be found so easily. Another method is to calculate the normals based on looking at points to the left and right of the point in question. Then find a line that goes through those two points and use a line perpendicular to it that goes through the center point to be the normal direction. This method would be difficult in the point set implemented thus far because the points resulting from the itkBinaryMaskToNarrowBandPointSetFilter are not uniformly spaced around the border. Because of the way the points are placed, in some cases two points on the same normal vector, this method would only work if the point set were constructed more carefully. Reducing the width of the filter causes there to be empty areas in the border that would not work effectively. Creating the ellipse by using equations to map out the desired region seems like it could also be one solution.

Aside from the creation of the point set normals, their incorporation into the registration will require either adding them to a point set class as an inherited class, or having them as a separate object. There will also need to be an efficient way to transform the normals that could work interchangeably with the various transforms. Those tools could be a part of the new point set class.

I was unable to get the program that the meshes are fed into to calculate the normals to work on my computer. Errors seem to be stemming from the inclusion of "argtable" into the program. I tried installing an rpm of the program from the Internet but get the error "Missing Dependency: linux-gate.so.1 is needed by package argtable2". However, "linux-gate.so.1" is supposedly a virtual library and cannot be installed. It may have something to do with where the install is looking for it. In my case since I am calculating the normals this issue is not as critical.