[vtkusers] Extracting largest "image" object, filling, adding.

gte631d at mail.gatech.edu gte631d at mail.gatech.edu
Tue Jul 15 21:12:02 EDT 2003


Hi All,
Purpose: Extract the largest connected 3D object (from vtkImageData)

** I have first explained my problems.  The specific PYTHON code (for vtk4.0 is 
appended) 

Intermediate Problem 1: vtkImageSeedConnectivity - only searches in the +X 
direction.  The documentation states "If a seed supplied by the user does not 
have pixel value "InputTrueValue", then the image is scanned +x, +y, +z until a 
pixel is encountered with value "InputTrueValue". This new pixel is used as the 
seed ."  As far as I can tell, this is true only for +X.  The only way for me 
to get it to work is by looping in Y and Z setting seeds one by one in the 
following manner.
	seed_im.RemoveAllSeeds()
	seed_im.AddSeed(0, seed_Y, seed_Z)

Intermediate Problem2: vtkImageMathematics - Simply adding two binary images 
within a loop seems to be quesionable.  For vtk 4.2 the following will work
	math_im = vtkImageMathematics()
	math_im.SetOperationToAdd()
	math_im.SetInput1(canvasUCimg) # Needs to be set each time
	math_im.SetInput2(clusterimg) 
	math_im.Update() #!!! Will UPDATE AGAIN at seed_im.SetInput if "=" used
                          # => essentially like addind clusterimg TWICE!
	canvasUCimg.DeepCopy(math_im.GetOutput())

I would think that the above is working properly if in place of DeepCopy I 
could use the following
	canvasUCimg = math_im.GetOutput() 
but if I do, it looks like Update occurs twice. ie. Looks like 1+2 = 5 (1+2+2) 
once canvasUCimg has been through vtkImageSeed connectivity inside (back at the 
top of) the same loop.

** for vtk4.0 this is even worse.  I just use the following to manually add 
(very slow)
	CanvScl = canvasUCimg.GetPointData().GetScalars()
	ClstScl = clusterimg.GetPointData().GetScalars()
	for id in range (0, img_size):
		CanvScl.SetValue(id, CanvScl.GetValue(id)+ ClstScl.GetValue(id))

#-------------- Here is My Function Code -------------------------#

def vtkImageLargestObject(binimg, connect_val) :
    #1. Cast input image to Unsigned Char; invert if looking for
    #   largest ZERO object
    threshcanv = vtkImageThreshold() # will be recycled
    threshcanv.SetInput(binimg)
    threshcanv.ThresholdBetween(connect_val, connect_val)
    threshcanv.SetInValue(1) # retained through recycling
    threshcanv.SetOutValue(0)# retained through recycling
    threshcanv.SetOutputScalarTypeToUnsignedChar()
    threshcanv.Update()
    canvasUCimg = threshcanv.GetOutput()  

    #2. Next Initialize the Necessary intermediate objects
    seed_im = vtkImageSeedConnectivity()
    seed_im.SetInput(canvasUCimg)
    seed_im.SetInputConnectValue(1)
    seed_im.SetOutputUnconnectedValue(0)
    seed_Y, seed_Z = 0, 0
    seed_im.AddSeed(0, seed_Y, seed_Z)
    
    labelledimg = vtkImageData()#DeepCopy canvas with labelled clusters to this

    hist_im = vtkImageAccumulate()
    hist_im.SetComponentSpacing(1,0,0)# 1D histogram

    #3. Initialize the output objects
    lgstclustimg = vtkImageData() #DeepCopy temp output to this
                                  #because we recycle threshcanv

    #4. Initialize attributes/flags
    max_count = 0
    max_bin = 2 # The first found cluster (if exists) will be labelled 2
                # if no cluster is found then the output will be all-zero.
    filling_complete = 0
    img_extent = canvasUCimg.GetExtent()
    img_size = canvasUCimg.GetPointData().GetScalars().GetSize()

    #5. Now begin the cluster filling and bin counting (histogram) prodedure
    while 1: #outer loop for the event that the considered image has more than 
254 clusters
        fill_val = 1 # initialize the label val for each new (semi) canvas

        while 1: #perform cluster labelling as long as necessary
            seed_im.SetOutputConnectedValue(fill_val)    
            seed_im.Update()
            clusterimg = seed_im.GetOutput()
            cluster_rng = clusterimg.GetScalarRange()

            while cluster_rng[1] == 0 : # Update seed until a new cluster is 
found
                seed_im.RemoveAllSeeds()
                if seed_Y < img_extent[3] :# or extent
                    seed_Y += 1
                elif seed_Z < img_extent[5]:
                    #seed_Y has reached its extent for a given seed_Z. Start on 
a new Z-slice
                    seed_Y = 0
                    seed_Z += 1
                else : # All necessary seeds have been exhausted.
                    filling_complete = 1 # set the flag
                    break 
                                              
                seed_im.AddSeed(0, seed_Y, seed_Z)
                seed_im.Update()
                clusterimg = seed_im.GetOutput()
                cluster_rng = clusterimg.GetScalarRange()
                if cluster_rng[1] > 0 : # a new cluster has been found
                    break
            #END while (cluster_rng[1] == 0)

            # add the newly labelled cluster to canvasUCimg
            # vtkImageMathematics(), should be used only after its Update() 
method is fixed
            CanvScl = canvasUCimg.GetPointData().GetScalars()
            ClstScl = clusterimg.GetPointData().GetScalars()
            for id in range (0, img_size):
                CanvScl.SetValue(id, CanvScl.GetValue(id)+ ClstScl.GetValue(id))
                             
            # Update label value for next cluster
            fill_val += 1

            if (fill_val >= 255) | filling_complete: # max possible fill_val is 
254
                print "Filling Rotation Complete!"
                labelledimg.DeepCopy(canvasUCimg)
                break
        #END while (1): inner while loop
        # While loop could have been exited by one of two ways:
        #   A. Finished labelling all clusters (and there were <= 254 clusters)
        #   B. Ran out of labels
        PrintImageScalars(labelledimg)

        # 6. Histogram of the labelled output
        labelrng = labelledimg.GetScalarRange()
        NumPixelLabels = labelrng[1] + 1 # REMEMBER: this function is 
written/optimized
                                         # for vtkImageData that is obtainted 
from
                                         # vtkImageSeedConnectivety()
                                         # because it has UnsignedChar scalars 
the min will be 0
                                         
        hist_im.SetInput(labelledimg)
        hist_im.SetComponentOrigin(2,0,0)
        hist_im.SetComponentExtent(0,labelrng[1]-2,0,0,0,0)
                        #bin each pixel-label seperately we don't want and
        hist_im.Update()
        histimg = hist_im.GetOutput()
        # we now have a 1-D array of hist bin counts.

        # 7. Get Max count and compute an output candidate
        histrng = histimg.GetScalarRange()
        if histrng[1] > max_count: #proceed and update maxcount
            max_count = histrng[1]
            histdims = histimg.GetDimensions()
            for scalar_id in range (0, histdims[0]):
                if (histimg.GetPointData().GetScalars().GetValue(scalar_id) == 
max_count):
                    max_bin = scalar_id + 2
                    break
           # At this point we have the (first) bin value of the largest cluster

            threshcanv.SetInput (labelledimg)
            threshcanv.ThresholdBetween(max_bin,max_bin)
            if connect_val == 0 :
                threshcanv.SetInValue(0)
                threshcanv.SetOutValue(1)
            threshcanv.SetOutputScalarTypeToShort()
            threshcanv.Update()
            lgstclustimg.DeepCopy(threshcanv.GetOutput())# Use deepcopy
                      # because threshcanv may be recycled below.
      
        if filling_complete : #exit the outer while loop
            print "Largest Object (if any) has been found!"
            return lgstclustimg

        #if still not broken out: 
                 #we have another round of filling/labelling to do
        #ie. There "may" still '1' voxels available

        # 8. Refresh the canvas with only unlabled 1 voxels    
        threshcanv.SetInput (labelledimg)
        threshcanv.ThresholdBetween(1,1)
        if connect_val == 0 :
            threshcanv.SetInValue(1) # reset because could have been changed
            threshcanv.SetOutValue(0)# if connect_val == 0.  Else: redundant
        threshcanv.SetOutputScalarTypeToUnsignedChar()
        threshcanv.Update()
        canvasUCimg = threshcanv.GetOutput()
    #END while(1): outer while loop
#---------------- End Function Code ---------------------#

I need help with **Adding** the images efficiently.

Any ideas suggestions will be appreciated.

Thanks,
Sarju


-- 





More information about the vtkusers mailing list