-1

Now there is a 3D volume stored in 3d array, a point is given in the 3d volume. How to calculate the minimum cross section of the volume through the point?

Silver0427
  • 63
  • 5
  • 1
    what is the size of the dataset? could easily do a numeric solve if the data is reasonably sized or processing time is not critical – Sneaky Polar Bear May 30 '18 at 12:34
  • @SneakyPolarBear The size of dataset is 180*512*512 3d array, background is 0 and volume is 1 and the processing time is not critical. Could you please show me how to do it? – Silver0427 May 31 '18 at 07:25

1 Answers1

2

you may be able to do some fancy footwork with a distance transform, but not sure how I would do that atm, so just gonna state it and leave it as I said I would give a numerical methodology

Also, be careful saying "processing time is not critical" as everything gets relative when you are doing a numerical solve on a data set that big. IE small changes in your inner loop could change processing time from minutes to hours to days...

Everything I am describing below is "off grid", what I mean by that is: you have a structured point cloud, but the spatial vector generation and propagation is in free space. So in order to, for instance, test a point for inside or outside along a given free space vector, you will look up the nearest voxel coordinate inside the structured cloud. Usually, the easiest way to do this is simply to scale things to "floor" to the nearest resolution unit. Although this is easier if you use integer resolution and can actually use floor, you can also just do it with division, rounding and re-multiplication along primary axes.

First step would be to set up a function to generate an arbitrary cross section. To do this, you will take a vector (A), then generate a perpendicular vector (C) to this vector (A) (this can be done in a number of ways, but easiest is to just create another non-parallel vector (B) and cross that with the original to create a perpendicular vector (C). For vector (B) you could just use (0,0,1) or (0,1,0) *using the former unless the test vector is parallel to (0,0,1) in which case you use the latter. Now that you have the vector (A), and a vector (C) in its perpendicular plane, you need to test your data along that vector (C). Basically, test origin point and then move along the vector (C) by your smallest resolution unit and keep testing until you reach the edge of your dataset or if you reach the edge of your object depending if you have a concave hull or not. Then rotate your vector (C) about vector (A) and radially sample the entire plane. That rotation step size should probably relate the the minimum arclength of sample data furthest from your sample point (you can probably just set some arbitrarily small radial sample distance, but this is center of processing loop, so stuff is most expensive at this step). It is important that you use some method to normalize your area after this step, as you will be repeat counting a lot of area. There are a number of ways to get the inside area of an arbitrary closed shape. The one that comes to mind is to just calculate triangle area from the last "in point" for each test vector (ie the edge) to the next and sum all those triangles to get the final area. If you are doing complex (non-concave shapes), you would want to get in and out edges and do trapazoidal area calculations for second and third edges.

That function will be run on a series of seed vectors that you generate at the beginning and emanate from the test point. These seeds will be a radial distribution of vectors to make a pretty reasonable sphere (depends on the object you are looking at but basically the seeding step is to ensure you avoid most large local minimums (so you may be able to have higher or lower seeding resolution depending on your object).

After you test all your seeds, for each seed, get its cross section, then attempt small rotations on its two perpendiculars and do gradient descent to find a local minimum. After testing all seeds, take the best solution.

Sneaky Polar Bear
  • 1,611
  • 2
  • 17
  • 29
  • Thanks for the detailed answer. I'll try the Least-squares – Silver0427 Jun 01 '18 at 06:19
  • did u mean I have to use a function to return the area of the cross section then use optimize operations to get the minimum value? – Silver0427 Jun 01 '18 at 08:45
  • yes, the initial seeding is just arbitrary guesses. An optimization (gradient decent or otherwise) function will repeatedly call the cross section area function while making small changes to that test vectors orientation. The easiest numerical way to do gradient decent is to just rotate a little one way, if it gets better, do it again, if it gets worse, rotate half in other direction. keep doing that until u are below a min step size (then do the other axis). There are a lot of better styles of this tho – Sneaky Polar Bear Jun 01 '18 at 12:10
  • 1
    Actually the 3D data displays the shape of a segment of aorta, so it's irregular. I think it's a little bit hard to present the area of cross section. In this case, what should I try? – Silver0427 Jun 02 '18 at 06:49
  • That actually makes it much easier, as you can limit your max projection distance to the reasonably largest diameter that the aorta could be which i assume is much smaller than the full dicom volume – Sneaky Polar Bear Jun 02 '18 at 15:56