at_plane (3D equivalent of at_line)

I was looking into ways to create a 2D image from an arbitrary cut on a 3D image stack formed by axes aligned slices. This would be equivalent to -at_line which returns a line in a (2D) image. Is this available on gmic, maybe under another name? I couldn’t find it.

Given a 3D point and a normal at it, that defines an oriented unbounded plane at which the closest points in the cut and on the positive side of the plane will be projected into thus creating a 2D image (with possibly ‘empty’ out of bounds parts). How could this be done with the current operators in gmic, if not there yet?

-Alex

That’s an interesting command to add indeed. If I have some time today, I’ll try to find out how I can do this.

Many thanks for the interest and willingness to implement it. It can be pretty useful :slight_smile:

No doubt you’ll come up with the best solution, but one that occurs to me is to rotate a 2d image (warp plane) holding x,y,z points then simply -warp using it. That way you can have nice interp. A moving origin could be added/subtracted before/after rotation.

For example, I have a quaternion matrix based point rotation to rotate x,y,z points around an arbitrary axis (angle is in range 0-1, could easily be converted to radian):

gcd_rotate : -skip ${1=0},"${2=1}" # $1 = angle, $2 can be a 3 channel vector
  -repeat $! -l[-1] -i ${2}
    -if {s!=3} -k[0] x={3^-0.5} y=$x z=$x
    -else --norm. -/[-2,-1] x={@0} y={@1} z={@2} -rm. -endif
    a={$1*pi} C={cos($a)} S={sin($a)} O={1-$C}
    -mix_channels (\
    {[$C+$x*$x*$O,$x*$y*$O-$z*$S,$x*$z*$O+$y*$S]};\
    {[$y*$x*$O+$z*$S,$C+$y*$y*$O,$y*$z*$O-$x*$S]};\
    {[$z*$x*$O-$y*$S,$z*$y*$O+$x*$S,$C+$z*$z*$O]})
  -endl -done

Then you can do something like (forgive the hastily written code):

100%,100%,1,1,"x-w/2"
100%,100%,1,1,"y-h/2"
100%,100%,1,1,0
-a[-3--1] c
-gcd_rotate. -0.3,(0^0^1)
-sh. 0 -+. {0,w/2} -rm.
-sh. 1 -+. {0,h/2} -rm.
-warp.. .,0,2,0
--k..
-c 0,255

The (0^0^1) being the axis.

Edit: if this is useful I could write a proper solution.

I’ve since discovered the -rotation3d command which could do the same job!

Thanks for sharing the code. I am not sure though how rotate3D can give us a 2D image which is an arbitrary cut on a 3D image stack. Sorry, it is probably my non-understanding of the steps in the code.

An illustration movie of possible 2D cuts on a cube follows (the image stacks I have in mind are not cubes but the cut possibilities are the same, I think):

-rotation3d provides a 3x3 rotation matrix, which can be used to rotate x,y,z points of a 2D “warp” field. It’s a shorter way of doing the above… sadly I probably won’t have time to explain/demonstrate in full until the weekend. If you can wait until then, or somebody else has a complete solution before me then all is good :slight_smile:
It can absolutely be done, the only question is which way… I assume a volumetric input image?

Edit: perhaps I’ve read what you want wrong - is it simply a 2D ‘slice’ you want?

That’s correct, it is simply a 2D slice that I want. A slice that cuts a given volumetric image. Pixels in the slice would have the intensities of those voxels in the volumetric image intersected by the slice plane.

I a straightforward way, if Ax + By + Cz + D = 0 is the plane equation describing the slice, one look for those discrete grid points (u,v,w) comprising the volumetric image that are ‘on’ the plane or ‘almost on’ the plane, i.e. points for which Au + Bv + Cw + D ~= 0 (or < epsilon = sqrt(3)/2 - half unit cube diagonal ). The intensity (color) of those points define the intensity of pixels on the plane. One can then do an orthogonal projection of the 3D plane to obtain a 2D image. But I tend to believe there might be better, faster ways to do this - the clipping plane in graphics might do this.

An equivalent in Matlab:

http://www.mathworks.com/help/matlab/ref/slice.html

@alosca Just curious: what would be your application for such a function?

@afre The application is to extract planes showing the interior parts of a 3D volumetric data. See the Matlab link I added in my previous post.

Is the pink 3dplane a solution?

Had some time today to think about it, and ended up with a new command -at_quadrangle in latest G’MIC 2.1.0.

   at_quadrangle:
                        x0[%],y0[%],x1[%],y1[%],x2[%],y2[%],x3[%],y3[%],_interpolation,
                          _boundary_conditions |
                        x0[%],y0[%],z0[%],x1[%],y1[%],z1[%],x2[%],y2[%],z2[%],x3[%],y3[%],z3[%],
                          _interpolation,_boundary_conditions

        Retrieve pixels of the selected images belonging to the specified 2d or 3d quadrangle.
        'interpolation' can be { 0=nearest-neighbor | 1=linear | 2=cubic }.
        'boundary_conditions' can be { 0=dirichlet | 1=neumann | 2=periodic | 3=mirror }.
        
        Example: [#1]  image.jpg params=5%,5%,95%,5%,60%,\
                   95%,40%,95% --at_quadrangle $params polygon.. 4,$params,0.5,255

This is a bit more flexible than the original request, as this command is able to ‘crop’ 2d or 3d image data from an arbitrary quadrangle.
For instance :

$ gmic volumes/reference.inr --at_quadrangle 0,0,0,100%,0,100%,100%,100%,100%,0,100%,0

gives these two images (input image first, then crop by quadrangle):

gmic_ref0

gmic_ref2

I’ve used one the diagonal plane to define the quadrangle.

I guess this should be easy to use this with custom coordinates to adapt it for a quadrangle that is a given plane with an orientation.

Also, I’ve added a new filter using command -at_quadrangle in the G’MIC plug-in. Does the same kind of crop, but limited to 2d images in this case.

What do you think ?

Thanks for working on this! I will download the latest and give it a try.

It can happen that the intersection of an arbitrary plane with the bounding box (a prism) of the volumetric image define only 3 or more than 4 points in the bounding box edges. Can I give the quadrangle coordinates outside the image?

Yes, you can specify quadrangle vertices outside the initial image domain, like in this example:

$ gmic -sp lena --at_quadrangle -50%,-30%,150%,-30%,40%,120%,60%,120%

That’s even the only case where the parameter boundary_conditions has a meaning!