Filling a closed curve ROI from coordinates of the vertices

Hello,

I have a series of coordinates which define some points that are the vertices of a closed curve.
Using 3dUndump, I am able to make a 3D dataset where those points have value=1 and the rest is 0 (see attached image, where the points are marked in red).

Is there an AFNI program to assign 1 to all the voxels inside the closed curve defined by those vertices? (Similar to what could be done with “Flood” on the drawing plugin)

Thank you,
Giuseppe

Screen Shot 2022-11-21 at 8.12.21 PM.png

I don’t think a command line program could do this generally.

As you mention, though, this is doable in the plugin.

I’m curious where the points came from? Is there some other step in their creation that could be used to help this goal, perhaps?

–pt

Thank you for your response.
I tried 3dRowFillin on different directions, but the result is not optimal.

The points were manually determined using another software, ImageJ.

Giuseppe

Well, OK. If your points are all in a single plane, as in that image, here is a script to connect successive points (defined as a 3-column text file of XYZ locations, and the first line should be repeated as the last, so that the shape gets closed.

I used the following file as the “points.txt” below:


-6.2    6.2   43.8
    6.2    6.2   43.8
   16.2    6.2   41.2
   23.8    6.2   36.2
   28.8    6.2   26.2
   21.2    6.2   18.8
   11.2    6.2   21.2
   -8.8    6.2   23.8
   -3.8    6.2   11.2
    6.2    6.2    6.2
   -3.8    6.2   -1.2
  -16.2    6.2    3.8
  -28.8    6.2    8.8
  -31.2    6.2   21.2
  -28.8    6.2   38.8
  -13.8    6.2   41.2
   -6.2    6.2   43.8

… and the input dataset (for setting the grid), was the “mask_epi_extents+tlrc.HEAD” in the output FT.results directory of the AFNI Bootcamp s05* output script. but it should probably work in many datasets with a reaonable field of view, just as an example.

The images show the points-to-be-connected as white underlay points, and then lines that are created to connect them (just looping over all successive pairs of points in the file, walking with a reasonably small step size (based on min voxel dim) to the next point), and then the planar infill of that.

Hope that might be useful. This would not extend to a case of arbitrary points in 3D—it depends on them being in a single plane, so that a closed curve can be drawn.

–pt


#!/bin/tcsh

# this program tries to connect points with a line.  The file "

set dset_grid = mask_epi_extents+tlrc.HEAD      # dset with grid to use
set file_pts  = points.txt                           # input text: XYZ coord
set dset_lines = dset_connect_the_dots.nii.gz        # output dset

# ----------------------------------------------------------------------

set npts   = `cat "${file_pts}" | wc -l`       # number of points in file

# set step size from voxel dims
set ad3    = `3dinfo -ad3 "${dset_grid}"`      # input voxel dims
set minad3 = ${ad3[1]}
foreach ii ( `seq 2 1 $#ad3` )
    if ( `echo "${ad3[$ii]} < ${minad3}" | bc` ) then
        set minad3 = ${ad3[$ii]}
    endif
end
set step   = `echo "scale=5; ${minad3}/3.0" | bc`
echo "++ The step is: ${step} (from vox: ${ad3})"


# create file: temporary list of points for this line
set file_tmp = _tmp_file_line.txt
set dset_tmp = _tmp_file_line.nii.gz
printf "" > ${file_tmp}

foreach nn ( `seq 2 1 ${npts}` )
    
    @ mm = ${nn} - 1

    # XYZ coords of points
    set xyz0 = `sed -n ${mm}p "${file_pts}"`        # start point
    set xyz1 = `sed -n ${nn}p "${file_pts}"`        # end point

    # dist between start and end points
    set dist = `echo "scale=5; sqrt(($xyz1[1] - $xyz0[1])^2 + ($xyz1[2] - $xyz0[2])^2 + ($xyz1[3] - $xyz0[3])^2)" | bc`

    # how many steps to take; shouldn't need to bother about
    # remainder, bc both endpoints are already included
    set nstep = `echo "scale=0; ${dist}/${step}" | bc`
    
    echo "++ The dist and nstep are:  ${dist} and ${nstep}"

    set coord = ( $xyz0 )
    foreach ss ( `seq 1 1 ${nstep}` )
        foreach ii ( `seq 1 1 $#ad3` )
            set coord[$ii] = `echo "scale=5; ${coord[$ii]} + ${step}*($xyz1[$ii] - $xyz0[$ii])/${dist}" | bc`
            # start walkin'
            printf "%10.6f %10.6f %10.6f\n" ${coord[1]} ${coord[2]} ${coord[3]} >> ${file_tmp}
        end
    end


    echo "-----------------"
    cat  "${file_tmp}"
    echo "-----------------"

end

# remove repeated locations
cat  "${file_tmp}" | uniq > FINAL.txt

#create a starter file with the given points
3dUndump                                                                \
    -overwrite                                                          \
    -xyz                                                                \
    -orient  RAI                                                        \
    -master  "${dset_grid}"                                             \
    -prefix  "${dset_tmp}"                                              \
    -datum   short                                                      \
    FINAL.txt

# fill holes in any plane (I think)
3dmask_tool                                                                  \
    -fill_dirs   AP                                                          \
    -fill_dirs   IS                                                          \
    -fill_dirs   LR                                                          \
    -fill_holes                                                              \
    -overwrite                                                               \
    -inputs      "${dset_tmp}"                                               \
    -prefix      "${dset_lines}"

exit 0

pts_with_lines.jpg

pts_filled.jpg

… and I notice now that the infilling did a bit more than desired, based on relatively convexity… Have to ponder about that.

–pt

Thank you so much for working on this, Paul, I truly appreciate your support!
I was able to run your code and get some results. I noted that the final mask output varies based on the order of the coordinates.
Specifically, I tested it on a case where my points were laying on a Z plane (so the Z coordinate is fixed for all the points).
I tried three different scenarios:

  1. CASE 1: the coordinates in the input file were exactly as I was given (i.e., in the order that the drawer traced them on ImageJ, in this case it was clockwise. In my project, the tracer will trace the coordinates either clockwise or counterclockwise). I just copied the coordinates of the first point at the end of the file as you suggested for your code.

292 297 80
294 298 80
296 297 80
298 294 80
297 292 80
297 292 80
294 292 80
292 292 80
292 294 80
292 297 80

  1. CASE 2: the X coordinates were sorted from low to high, and for each X, the Y coordinates were sorted from low to high.

292 297 80
292 292 80
292 294 80
294 292 80
294 298 80
297 292 80
297 292 80
296 297 80
298 294 80
292 297 80

  1. CASE 3: the Y coordinates were sorted from low to high, and for each Y, the X coordinates were sorted from low to high.

292 292 80
294 292 80
297 292 80
297 292 80
292 294 80
298 294 80
292 297 80
296 297 80
294 298 80
292 292 80

The results from the 3 cases are attached. To “fill” a convex polygon your program should have the X coordinates sorted from low to high, and, for each X coordinate, the Y coordinates should also be sorted.
The result from case n.2 was indeed the one that fitted my goal the best, although not exactly. As you can see in the attached file “goal”, I would like to have segmented also the voxels reported in green.
Anyway, thank you so much for making this program. I will try to work on it to achieve my final goal. Thank you for your help on this, really appreciate it.
Please let me know if you further update your program or if you have any suggestion. Thank you!

Giuseppe

Hi, Giuseppe-

Indeed, the order of points in the text file matters—they would have to be given in the order you want the line segments to be drawn.

To be clear, this is a veeery primitive algorithm for doing this… Maybe there is a cleverer way, the way that shapes could be made is so general, it is pretty tricky.

I think using the draw dataset plugin might be the best way to accomplish this, actually, even though it isn’t scriptable.

How many datasets do you have to process in this way?

-pt

I have ~5000 datasets, that’s why I am trying to use a script.

I just found out that in order to have the shape got closed by your program, the first line should be repeated as the last not once, but twice.
By doing this, the shape gets closed and the mask that I obtain matches exactly what I need. I used the original list of coordinates, which are reported following either a clockwise or an anti-clockwise direction around the polygon.

Thank you SO much! Maybe it would be useful to implement this program in AFNI, as other people might need to do a similar task, maybe!

Giuseppe

Hi, Giuseppe-

OH. 5000 datasets probably is more instances than you want to use in the GUI…

Glad that was useful for you.

Will ponder making that into a slightly more script/commandline form, indeed, then.

–pt

Paul’s solutions looks pretty good to me, but for 5000 ROIs, it might be a little slow, but that might not be too bad. We have something a little like this inside our Draw Dataset plugin - line propagation and then ROI filling, but it’s not exposed in a separate function. For reference, I think the line propagation is done in AFNI_3d_linefill() in afni_receive.c and DRAW_2dfiller in plug_drawdset.c. The filling is done in 2D in the slice orientation in which it was drawn. Unless the ROI is blobby, filling in all 3 directions, as noted in Paul’s script, may not be what you want, so you may want to tweak that as an option for a particular direction.

In any case,ImageJ and its branch, Fiji, have plugins themselves for doing this kind of job. I’m not sure this will be better than the AFNI way, but you might consider that in the future.

https://forum.image.sc/t/segmentation-masks-from-rois-plugins/55582
https://imagej.net/plugins/masks-from-rois

Also, just FYI, there are the poly2mask and regionfill functions in Matlab for doing similar processing in 2D.
https://www.mathworks.com/help/images/ref/poly2mask.html