This page assumes you have obtained a copy of the 3DMA-Fiber code, compiled it and need to start using it. If you want to obtain a copy of the code, email Prof. Brent Lindquist at lindquis@ams.sunysb.edu.

The historical reference, 3DMA General Users Manual, gives much of the basic low level information on internal file types, etc. but has not been updated for new features and directions. This web-based primer exists to provide up-to-date information to orient the user to, and allow the user to run, the current version of 3DMA-fiber. This document consists of the following major sections. It relies on an example data set to illustrate the software capability.

• General orientation of workflow
• Directory set-up
• Download example data set
• Executing the code
• List of input files


• Operations performed on tomographic data
• Case 1.1: Visualize the tomographic data
• Case 1.2: Resize the tomographic data
• Case 1.3: Rotate the tomographic data
• Case 1.4: Image intensity histogram, porosity-threshold dependence
• Case 1.6: Cylindrical RoI construction
  Automatic construction of the perimeter of a circular region of interest (RoI) to impose on each slice in the data stack. Subsequent analysis is restricted to the interior of the RoI (fiducial polygon).
• Case 1.7: Fiber/void space segmentation
  
  • Case 1.7smpl: Simple segmentation
  • Case 1.7krig: Indicator kriging


• Operations performed on segmented data
• Case2.1: Visualize segmented files
• Post-segmentation clean-up procedures
  • Case2.3: Disconnected volume distribution
  • Case2.4: Disconnected volume cleanup
• Case 2.5: Medial axis construction
• Case 2.6: Surface area computation


• Operations performed on medial axis data
• Case 4.1: Visualize the medial axis
• Case 4.4: Fiber thickness distribution
• Case 4.10: Medial axis modification (trimming)
• Case 4.12: Fiber reconstruction and analysis


• Operations performed on simulated data
• Case 0.1: Generate a fan of simulated fibers for producing a length correction table
• Case 0.2: Generate a simulated trellis-fan test data set
• Case 0.3: Convert simulated data to "tomographic" (digitized, grey-scale) data.


• Operations performed on void space through fiber mat
• Directory set-up, download of example data
• Case 2.5v: medial axis computation
• Case 4.10v: medial axis trimming
• Case 4.1v: medial axis visualization
• Case 4.6vx: tortuosity x-direction
• Case 4.6vy: tortuosity y-direction
• Case 4.6vz: tortuosity z-direction
• References

Fiber analysis consists of 3 general steps:
   user provided greyscale, digitized data (assumed to be tomographic in origin - though any greyscale digitized image is ok)
   - is converted (segmentation) to binary, two phase (fiber/non-fiber) data
   - from which the medial axis is extracted.
   - Analyses performed on the medial axis extracts the fiber information.

These three general steps are accompanied by a fourth step to produce a fiber length correction table. This step may only have to be done once to provide a correction table appropriate for analysis of a number of fiber data sets. It may however have to be redone if there are significant differences in length, thickness or digitized resolution of the fibers.

Fig. 1:   Workflow overview for fiber analysis.

Fig. 1 provides the workflow overview for fiber analysis. There are, generally two parallel sequences involved. The first is an analysis of the user data, consisting of the general data sequence:

tomographic -> segmented -> medial axis.

simulated -> digitized (tomographic) -> segmented -> medial axis -> length correction table

Thus (in either streams) the analysis can be viewed as flow through the short sequence of data sets:

tomographic -> segmented -> medial axis

Notes:

• The simulated data set, which consists of a 3D angular fan of non-touching fibers of uniform length and thickness, is digitized and run through the same tomographic -> segmented -> medial axis analysis in order to accurately capture all algorithmic systematics that affect the user data to be ultimately corrected.
• The tomographic -> segmented -> medial axis flow sequence is a general sequence upon which all of our 3DMA software is based. The software algorithms are organized to reflect this sequence and the structure of this primer reflects this organization.

A number of operations are performed on each data type, in order to acquire parametric (decision) values needed for production of the next data type in the sequence. Thus 3DMA-Fiber consists of:

0. Operations performed to produce/operate on simulated data
1. Operations performed on tomographic data
2. Operations performed on segmented data
4. Operations performed on medial axis data

Historical, we refer to the operations as "cases"; each operation has a short-hand name consisting of the word "case" followed by 2 numbers ("X.Y"). The first digit "X" refers to the type of data being operated on, the second digit "Y" refers to the operation performed. Figures 2a-d graphically summarize the operations typically performed on each type of data set in fiber analysis.

Fig. 2a:   Summary of the operations, with their associate "case" labels, commonly used to produce/operate on simulated fiber data.	Fig. 2b:   Summary of the operations, with their associate "case" labels, commonly used to operate on the tomographic data in fiber analysis.

Fig. 2c:   Summary of the operations, with their associate "case" labels, commonly used to operate on the segmented data in fiber analysis.	Fig. 2d:   Summary of the operations, with their associate "case" labels, commonly used to operate on the medial axis data in fiber analysis.

• Operations acting on/producing simulated data are referred to as "case0" operations; on tomoraphic data are referred to as "case1" operations; on segmented data as "case2"; on medial axis data as "case4". (There are "case3" operations - but they are seldom used.)
• In each set of operations, there is one "case" that computes the next data type in the sequence; case1.7 algorithms convert tomographic data to segmented data; case2.5 algorithms convert segmented data to medial axis data.
• Case4.12 consists of the algorithms performing the final fiber identification and analyses.
• There are two separate medial axis data sets (Fig 2c). The first consists of two separate file types: distance information ("burn") files and medial axis files consisting of lists of medial axis voxels. The second medial axis file type consists of (three) files that organize the medial axis into a node/path network which encorporates both burn and medial axis information as well as elucidating the network structure of the medial axis.
• Complete information on the structure of the segmented and medial axis data files can be found in the appendix of the 3DMA General Users Manual.

- Create a directory for each fiber data set to be analysed. For this primer, we will use an example directory /nfs/3dma_fiber/Test/

  cd /nfs/3dma_fiber
  mkdir Test

We have developed a working subdirectory structure to organize the data files created by 3DMA-Fiber. None of this directory structure is hard-wired into the code, so you can adopt your own conventions. However, as some of the 3D intermediate files (the "burn" and voxel list "medial axis" files) are printed slicewise, the shear number of files produced can get large and we have found the follwing structure to work well.

- Organize the following subdirectories under /nfs/3dma_fiber/Test/:

   cd Test
   mkdir tom tom_rs tom_rot cases burn seg seg_smpl c_seg ma ma_t plots results

The subdirectory names are highly suggestive of what goes in each, here is a brief description of what this test example stores in each subdirectory:

tom		raw, tomographic files
tom_rs		resized tomographic files
tom_rot		rotated tomographic files
cases		input, output files
burn		distance (erosion or "burn" number) files
seg		segmented files obtained by the kriging method
seg_smpl		segmented files obtained from simple global thresholding
c_seg		post-processed (cleaned) segmented files
ma		medial axis files
ma_t		post-processed (trimmed) medial axis files
results		graphical output on geometrical characterization
plots		3D plot files

Notes:

• In general, tom, tom_rs, and seg_rot subdirectories are not all needed, depending, of course, on whether resizing and rotation of the data set is desired. They are included in this example to illustrate the capabilities of 3DMA-Fiber.
• Likewise, generally seg and seg_smpl directories are not both needed. As indicated the seg directory holds segmented data files produced by the kriging method, seg_smpl holds segmented data files produced by global thresholding. Except for exploratory purposes, only one or the other segmentation method (preferably kriging) will be used and only one of the two directories are needed. Both are included in the exmple for instruction.

Use the following link to download an archived (fbr-tom.tar) data set of example tomographic input files into your /nfs/3dma_fiber/Test/tom/ directory.

This archived data set is 4.1 MBytes in size.

The data set contains 10 files labeled 'fbr.zzz.gz' where zzz ranges from 001 to 010, containing 10 slices (the z-direction) of a 512 x 512 x 10 volume of an unwoven fiber mat. Each slice is an array of binary short integers representing X-ray attenuation coefficients. The binary shorts are stored big-endian (Motorola processors: SUN/SGI/MacIntosh machines) and may have to be byte swapped upon read-in if they are to be processed on a little-endian (Intel processors: PC machines) computer. This can be done by choosing the appropriate menu option on input data file type.

- Extract the tomographic files.
  unzip fbr-tom.tar.gz
  tar -xvf fbr-tom.tar

3DMA-Fiber is a command-line based code, querying the user with a scrolled menu of options. To avoid the agony of typing in answers, each functional step of 3DMA-Fiber can read these answers from a pre-prepared input file with answers provided in the expected order. To facilitate the preparation of input files, 3DMA-Fiber echos every menu and its user response into an output file, which can used later as an input file.

We have adopted a convention for naming the input/output file pairs associated with each functional step. The input file is labelled caseN.M.in and its associated output file is caseN.M.out where N and M are the integer answers provided to the first two menus presented by 3DMA-Fiber.

The first menu is

	Input Data Options
		simulated fan data (0)
		tomographic data (1)
		segmented data (2)
		burn data (3)
		medial axis data (4)
	Enter choice:

and refers to the (intermediate) set of data files on which the desired analysis is to be performed.

The second menu, which also takes a numerical answer, depends on the first menu selection.

As stated, the 3DMA-Fiber executable is built to read user responses either from the keyboard or from an input file (caseN.M.in), and produce an output file (caseN.M.out) containing a record of the input as well as a text record of pertinent output statistics. The output file caseN.M.out also contains information on the run time and total memory required to perform the computation. This provides critical "scale-up" information for performing runs on larger data sets. In addition, depending on the case, 3DMA-Fiber outputs intermediate data files, as well as files containing requested graphs and images.

To run 3DMA-Fiber with user input from the keyboard:

	3dma_fiber_linux > caseM.N.out

Typing answers from the keyboard can be tedious. In addition, error recovery from incorrect input is poor, input must be entered all over again. Once prepared however, an output file can be converted to an input file, edited to change responses, and rerun making subsequent running relatively easy.

To reuse an output file:

	mv caseM.N.out caseM.N.in
	vi caseM.N.in                 /* edit input file, change answers */
	stripcomm_linux -b < caseM.N.in | 3dma_fiber_linux > caseM.N.out

The stripcomm facility is provided for the purpose of running from an input file. It scans the input file searching for the next occurrance of the string ": ". The remainder of the line is taken to be the next required answer and is passed to the 3DMA-Fiber executable.

As the stripcomm command line is rather lengthy, we provide the shell script runcase to shorten typing. Thus

	runcase linux caseM.N

The following input files will be used to illustrate 3DMA-Fiber analysis on the downloaded example data set. Clicking on each link, loads the file into the web browser. The browser's File->Save As function can then be used to save the input file into the "cases" subdirectory.

Input file		Function
case1.1.in case1.1_rs.in case1.1_rot.in case1.1_fp.in		slicewise rasterplotting of tomographic data
case1.2.in		resizing tomographic data to a smaller block
case1.3.in		rotating tomographic data
case1.6.in		preparation of a circular cross section, cylindrical fiducial polygon
case1.4.in		histogram of attenuation coefficients, dependency of porosity on simple thresholding cut-off
case1.7smpl.in		segmentation by simple thresholding
case1.7.in		segmentation by indicator kriging
case2.1gv.in		subvolume (Geomview format) plot of triangulated surface separating pore-grain surface
case2.3.in		volume information on disconnected voids and fiber segments
case2.4.in		disconnected volume removal operations on segmented data
case2.5.in		erosion computation of media axis (Lee, Kasyap, Chu [1994] algorithm)
case2.6.in		specific surface area computation
case4.1.in		(sub)volume plot (Geomview format) of medial axis
case4.4.in		fiber thickness estimation (histogram of medial axis burn data)
case4.10.in		modification (trimming) operations on medial axis
case4.12.in		fiber tracing algorithm
case2.5slc.in case4.1slc.in case4.10slc.in case4.1slc_t.in case4.4slc.in		preparation and histogram of medial axis burn data for a single slice (fiber width estimation)

Case 1.1 plots images of the input tomographic data.

Notice our specification in the input file (case1.1.in) of the data type as "swapped binary short (sbs)" in order to convert the big-endian files provided in the example data set to little endian. If you are using a big-endian machine, then specify the "binary short integer (bs)" data type.

Test Example:
	To obtain slice-by-slice rasterfiles (2-D plots of each slice in the
	stack which are viewable almost everywhere)

	- download the input file case1.1.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.1

	- This will produce the output file case1.1.out and rasterfiles
	fbr.001ras.gz to fbr.110.ras.gz in the directory /nfs/3dma_fiber/Test/tom
	These files can be viewed by most 2D image viewing packages, e.g. "xv"
	or "gimp" (the latter comes with a RedHat Linux distribution).

	% cd ../tom/
	% xv fbr.*.ras.gz

Fig. 3:   Rasterfile image of the test example fbr.005.ras.gz. The light phase represents void, the dark is fiber.

Case 1.1 provides an option to narrow the dynamic range by avoiding low and high ranges of attenuation coefficients that are only sparsely represented.
In the case1.1.in file, replace the answer

	Select method for determining MIN and MAX values for scaling data
		volume          (v)
		slice           (s)
		slice  adjusted (sa)
	Enter choice: v
	Set global min and max from files?
	(if the answer is not yes, you will be asked to enter the values)
	(y,n(dflt)): y

	Select method for determining MIN and MAX values for scaling data
		volume          (v)
		slice           (s)
		slice  adjusted (sa)
	Enter choice: sa
	Enter minimum counts to determine adjusted values at low and high end
	of histogram. Separate values for the minimum counts may be specified
	at each end.
		low  end value: 100
		high end value: 100

Fig. 4:   Rasterfile image of the test example fbr.005.ras.gz using adjustment of the grey scale to heighten the contrast between the two phases.

Case 1.2 allows the user to cut out a smaller sub-block of the tomographic image for study. This is generally useful for: speedup; avoidance of bad data regions; and avoidance of air regions.

Note that case1.2 outputs the resized files in binary integer format, (and NOT in the same format as the original tomographic data). By default, most case1.X and case2.X related computations done in 3DMA-Fiber utilize four-byte storage (integer or float) arrays.

Test Example:
	- an example input file case1.2.in can be downloaded into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.2

	- This will produce the output files fbr.001.gz through fbr.005.gz
	in the directory /nfs/3dma_fiber/Test/tom_rs. These 5 slices contain
	a 291 x 291 x 5 subset of the original data. The case1.1.in file can
	be modified as shown in case1.1_rs.in to produce rasterplots of this
	reduced block.

Fig. 5:   Rasterfile image of slice 3 of a reduced subvolume of the test example data set.

Fig. 5 shows the resulting rasterfile image of slice 3 in the reduced subvolume of the test example data set. Slice 3 in the reduced set corresponds to slice 5 of the original.

For the example analysis sequence discussed here, we will make no use of these resized data sets. However, we commonly utilize this resize case to reduce 1024^3 data sets to a manageable subvolume. Viewing the rasterplots from case1.1 can be helpful in guiding this trimming.

Case 1.3 allows the user to perform simple 90 degree rotations of the data set and re-output the data in slices along the new z-direction. This is merely a convenience step, often used to reorient data originally sliced in CD-transverse planes into slices in MD-CD planes.

Note that case1.3 outputs the rotated files in binary integer format. See the analagously placed comment in case1.2.

Test Example:
	- Download the example input file case1.3.in into
	the /nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.3

	- This will produce the output files fbr.001.gz through fbr.291.gz
	in the directory /nfs/3dma_fiber/Test/tom_rot. These 291 slices are
	the rotated 291 x 5 x 291 version of the resized data in tom_rs.
	The case1.1.in file can be modified as shown in case1.1_rot.in to
	produce rasterplots of this rotated block.

As the rotated slices (of size 291 x 5) in this example are very thin, it is rather unsatisfactory to show a plot of the resultant rotated slices.

For the example analysis sequence discussed here, we will make no use of these rotated data sets. As mentioned above, we commonly utilize this case to conveninently re-orient the data slices.

Perusal of Figs. 3 and 4 shows that the fiber mat image has cylindrical cross section; data beyond a certain radius is spurious and should be ignored. 3DMA-Fiber allows the user to impose a cylindrical RoI with a simply connected polygonal cross section upon the image. Only the region interior to the cylinder will be analysed. For simple polygonal cross sections, the user can enter (x,y) coordinates of the perimeter vertices. This polygon is applied to each slice of the image; thus the volume to be analyzed is "cylindrical". We refer to the selected volume as the "fiducial volume" or the "fiducial polygon". (See section 4.5 of the 3DMA General Users Manual for details on specifying the x,y coordinates of the polygon vertices.)

Case 1.6 provides the ability to generate a polygonal approximation to a circular perimeter by specifying the center and radius of the circle. The polygonal vertices are output into a file. Note, the graphics viewer "xv" provides a nice way (through middle mouse button click) to find the desired center and radius values needed for a slice.

Test example:
	- download the input file case1.6.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.6

	- This will produce the output file case1.6.out and the coordinate
	file fp_fbr in the cases directory.  The case1.1.in file can be
	modified as shown in case1.1_fp.in to utilize the information in
	the fiducial polygon file to produce rasterplots of the tomograhic
	data showing the excluded region as demonstrated in Fig. 6.

Fig. 6:   Rasterfile image of slice 5 showing the region exluded by the fiducial polygon.

In all subseqent cases discussed in this primer, the fiducial polygon will be explicitly used (or implicitly when based upon previous use) to delimit the region analyzed.

The first major algorithm in the 3DMA-Fiber analysis is segmentation. Producing a distribution of the X-ray attenuation coefficients in the image is a decision-making tool useful for setting segmentation parameter values. Case 1.4 will produce a histogram plot (jgraph format) of the distribution density as well as a cumulative probability plot. The cumulative plot is interpretable as the fractional volume occupied by all voxels whose intensity lies below any X-ray attenuation coefficient value, and thus directly provides a measure of porosity versus cut-off threshold for the global (simple) segmentation method.

Test example:
	- download the input file case1.4.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.4

	- This will produce the output file case1.4.out and the jgraph file
	  att.jgr in the directory /nfs/3dma_fiber/Test/results
	  The jgraph file can be processed and viewed as follows:

	% cd ../results/
	% jgraph < att.jgr > att.ps
	% ggv att.ps

The jgraph program is publically available software that is distributed with 3dma code. It produces postscript output which can be printed, included in documents, or viewed by a postscript interpreter (e.g. ghostscript). Fig. 7 shows the resulting plots obtained in the above example output file att.ps.

Fig. 7:   Histogram of X-ray attenuation coefficients (dashed line, right axis - not normalized to unity) and the cumulative probability (solid line, left axis) interpreted as void fraction (porosity).

As is typical with X-ray CT images of bi-phase materials having sufficient X-ray attenuation contrast, the plot is bimodal, the lower mode corresponding to the lower attenuating phase (void space) and the upper mode to the higher attenuating phase (fiber). The void peak corresponds to attenuation coefficient values near 0, the fiber peak is around 2000. The wide dynamic range of attenuation coeffient values, extending from -6,000 to 32,000, is responsible for the lack of void-fiber contrast seen in Fig. 3.

Simple (global) thresholding interprets each voxel having attenuation coefficient value below a threshold as phase "0" and every voxel above as phase "1". It is only reasonably useful as a fast segmentation in order to do a quick analysis. Time permitting and quality counting, the kriging based segmentation case (case1.7) should be used.

Choosing an example threshold value of 1,000 (based on the histogram of case1.4, Fig. 7) results in the following test example run.

Note that we have coded simple thresholding with two thresholds. Read the comments in the input files carefully. If the lower threshold is set below the minimum attenuation coefficient value of the data, the lower threshold is invalidated; similarly if the upper threshold exceeds maximum data values, it is invalidated. Attenuation coefficient values lying within the lower to upper threshold values are always segmented to unity; those outside to zero.

Test example:
	- Download the input file case1.7smpl.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.7smpl

	- This will produce the output file case1.7smpl.out,
	  intermediate (segmented) data files fbr.XXX.gz in the 
	  /nfs/3dma_fiber/Test/seg_smpl directory, and rasterfiles 
	  fbr.XXX.bwras.gz also in the ../seg_smpl/ directory.

	  The rasterfiles provide a graphical view of the segmentation 
	  results for each slice and just echo** the output of the
	  fbr.XXX.gz data files.

	% cd ../seg/
	% xv fbr_smpl.*.ras.gz

Fig. 8:   View of the test example simple segmentation results for slice 8 in the 3D image stack. Void space is white; fiber black.

As demonstrated in Fig. 8, simple segmentation is extremely prone to "speckled noise" and is best employed as a fast, rough, segmentation method for preliminary exploration of data.

The indicator kriging method requires input of two thresholds, T0 and T1. T0 (T1) establishes a lower (upper) limit below (above) which any voxel can reasonably be identified as phase 0 (1). All other voxels (with coefficient values in the range [T0,T1]) have their phase determined by indicator kriging. In the test example below, we use T0 = 500, T1 = 1500.

Test example:
	- download the input file case1.7.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case1.7

	- This will produce the output file case1.7.out
	  (in the cases directory), the intermediate data files
	  fbr.XXX.gz in the /nfs/3dma_fiber/Test/seg directory, and
	  kriging summary raster files fbr.XXX.ras.gz and fbr.XXX.bwras.gz
	  also in the seg directory.

	  The *.ras.gz rasterfiles provide a visual synopsis (Fig. 9) of
	  the kriging decisions for each slice. The *.bwras.gz rasterfiles
	  show only the final segmented images of each slice (i.e. show only
	  the bottom left plot in Fig. 9).

	% cd ../seg/
	% xv fbr.*.ras.gz

Fig. 9:   Kriging summary image for the 8'th slice in the test example.   (upper left) The tomographic data.   (upper right) White (black) voxels are identified as phase 0 (1) through the threshold T0 (T1). Yellow (red) voxels are identified as phase 0 (1) via kriging.   (lower left) The final segmented result.

The *.bwras files are plotted in this example only for completeness. The *.bwras files are unecessary if the *.ras files are kept. However the *.ras files require much more disk space to store, and most can be deleted. (One or two can be retained as demonstrations of the quality of the segmentation.) The *.ras files can be recovered from the seg/fbr.XXX.gz data files by running case 2.1, and specifying rasterfile slice output. The *.ras files cannot be recovered without re-running case 1.7.

With the image segmented, all further processing steps utilize the segmented data volume file. Complete information on the structure of these files is found in section 9.1 of the 3DMA General Users Manual.

No segmentation procedure is perfect - some level of segmentation errors result. In some respects, the most troubling of these are small isolated phase "blobs". Isolated blobs of phase 1** can be especially troubling if medial axis analysis of phase 0** is to be performed, as the isolated phase 1 blobs are cavities in phase 0 and will cause the medial axis of phase 0 to consist of a network of curved paths and closed surfaces. The closed surfaces surround the phase 1 blobs. A mixed path/surface medial axis network is difficult to analyse. The essential feature of case 2.4 is the ability to "delete" isolated blobs of a given phase (by converting them to the opposite phase). We use a volume threshold to decide what blobs to convert. In order to determine volume thresholds for isolated blob conversion of each phase, case2.3 should be run first. Case 2.3 provides both text printout and jgraph histogram of the distribution of isolated blob sizes. Unless there is some interest in retaining blobs above a certain size, in practice we convert all isolated fiber and all isolated void blobs. While isolated void spaces may be real in some fibers, their existence is usually related to errors. Floating (isolated) fiber "segments" are unphysical, except in one instance. The boundary of the imaged region may cut through a projecting fiber, leaving the fiber "tip" as a visible protrusion into the imaged volume, but losing the connectivity information of the fiber tip. In such cases, we consider it a small error to "delete" the protruding tip (by converting it to void phase).

Test example:
	- download the input files case2.3.in and case2.4.in
	into the /nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case2.3

	- This will produce the output file case2.3.out and a jgraph
	  file dv.jgr in the /nfs/3dma_fiber/Test/results directory.
	  The distribution ef volume sizes (in numbers of voxels) ordered by
	  volume for both disconnected pore and grain volumes is printed in
	  the case2.3.out file - as shown here.
	
	  To delete ALL blobs but the largest of each phase, volume thresholds
	  (e.g. 10 for phase 0 (air) and 1000 for phase 1 (fiber) can be
	  established simply by examining this printout. The chosen thresholds
	  can then be entered into case2.4.in.
	
	% vi case2.4.in
	% runcase linux case2.4

	- This stores a "cleaned" segmented volume file fbr.XXX.gz and raster
	  file plots fbr.XXX.bwras.gz of these in the /nfs/3dma_fiber/Test/c_seg/
	  directory.
	
	% cd ../c_seg/
	% xv fbr.*.bwras.gz

Fig. 10:   Slice 8 after disconnected volume clean up of the image.

Fig. 10 displays slice 8 after clean up. Although fibers may look disconnected in a single slice, after case 2.4 is run all fibers in the (fiducial polygon defined volume of the) 3D image are connected.

There are other correction options in case2.4.in. One is an attempt to correct arc (ring) artifacts that show up in the tomographic images as (arcs of) alternating black/white rings centered on the axis of rotation of the image and may result in (arcs of) rings of misidentified phase. This correction algorithm has had limited application. The remaining options allow for some type of fiber surface smoothing, and consist of two variations on "majority filtering" and the standard morphological closure operation. Unless there is specific reason to believe fiber surface roughness is affecting the analysis, it is probably better to leave these options non-functioning.

Analysis continues with the processed segmented files in the /nfs/3dma_fiber/Test/c_seg/ directory.

2D slice-wise rasterplots of segmented files have been encountered in the previous sections. Case 2.1 provides an option to produce rasterplots of (cleaned-up) segmented files by requesting the "Rasterplot (2D)" option in the "Image format options" sub-menu.

3D visualization routines have also been included to provide plots of the fiber interface (surface). The pore-grain surface can be plotted as either: a surface composed of (rectangular) voxel faces (producing a plot which is a 3D analog of the 2D rasterplots) in either Geomview or Inventor format; or as a triangulated surface obtained via the Marching Cubes algorithm [Lorensen and Cline, 1987]. This latter surface plot is available only in Geomview format and is illustrated in the test example below.

Test example:
	- Download the input file case2.1gv.in
	  into the /nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/ 
	% runcase linux case2.1gv

	- This will produce the output file case2.1gv.out and the Geomview
	  files fbr_surf.list and fbr_vol_axes.vect in plots/.

	% cd ../plots
	% geomview fbr_surf.list

The output from Geomview is shown in Fig. 11. Note we have only plotted a subvolume of the image. 3D plots can involve very large plot-files!

Fig. 11:  Geomview graphics view of fiber surface of the [200,299] x [200,299] x [1,10] voxel subvolume of the sample obtained by the Marching Cubes Algorithm

The processed segmented files are used to compute the medial axis of the fibers using the LKC [Lee, Kashyap, Chu; 1994] algorithm. In addition, an integer distance value (referred to as a "burn" or "erosion" number) is associated with each fiber voxel. This distance is the max norm distance from the voxel to the closest fiber surface voxel in units of voxels. (A fiber surface voxel has a "burn" number of 1.)

Test example:
	- download the input file case2.5.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case2.5

	- This will produce the output file case2.5.out,
	  intermediate burn number data files fbr.XXX.gz 
	  in the burn/ directory and intermediate medial 
	  axis data files fbr.XXX.gz in the ma/ directory.

The burn number files record distance information for all fiber voxels. All void voxels are assigned a burn number of 0. The medial axis data files record locations of the medial axis voxels only. Complete information on the structure of these files is found in section 9.3 of the 3DMA General Users Manual.

Case 2.6 computes the specific surface area (in units of voxel length - typically microns) of the fiber surface.

Test example:
	- download the input file case2.6.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case2.6

	- This will produce the output file case2.6.out in which the
	  specific surface area value is recorded.

	  e.g.
	  Volume-to-surface ratio 77.7266  Specific surface area 0.0128656

	  The units are length and length^{-1} respectively, where length
	  is in terms of the physical units specified in the input file.
	  In this example, the units of length are in microns.

Any medial axis algorithm is sensitive to surface noise resulting in dead-end paths that have no relevance to the geometrical features under investigation. These should be removed before the axis is utilized either to characterize the object or as a search path for feature location. Isolated pore spaces will produce isolated "pieces" of the medial axis which may or may not be of interest. Using terminology from graph-theory to refer to components of the media axis, case 4.10 provides the user with the ability to trim features from the medial axis. The major components of this trimming are:

• Remove isolated components corresponding to isolate pores: These come in two types, isolated paths, and single isolated medial voxels (paths of length 1 voxel). A path length threshold is used to determine what isolated paths to keep, a distance measure (burn number) threshold is used to isolated voxels.
• Remove dead end (branch-leaf) paths: Controls to determine a threshold to decide which to remove depend on path length and whether they exit through the boundary. For fiber studies, we recommend removing all branch-leaf paths whose length does not exceed the diameter (or perhaps twice the diameter) of a typical fiber. Branch-leaf paths exiting the boundary howewver, should probably not be removed.
• Needle eye-paths are paths which loop, i.e. return to the same branch point in the network from which they left. Unless it is known that there are closed loop paths in the pore space, all occurances of needle eyes are probably small path segments produced by local segmentation effects and should be removed.
• Reduce "surface remants": The medial axis consists of a network of paths that join at "branch points". In a digitized medium, a "branch point" may in fact not be a single voxel, but a set of voxels, which we refer to as a branch cluster. Ideally the set of voxels forming the branch cluster has genus 1 (one object, no holes or embedded cavities). In fact, it is possible to produce branch clusters that have holes. We refer to such clusters as "surface remnants", a phrase which is somewhat of a misnomer. 3DMA-Fiber can attempt to simplify the structure ot branch clusters, reducing them to genus 1 objects if possible (not guaranteed).
• Merge close clusters: one goal of the application of medial axis to studying fiber structure is to identify places where fibers bond. Unfortunately a bond site will probably involve at least two, and possibly more, closely connected branch clusters. The "merge close cluster" algorithm attempts to identify paths that link clusters belonging to a single fiber-fiber bond site by performing one of a number of test variants that compare the path length to local estimates of the distance to the fiber surface. If a path is judged to link two such close clusters, the path is "removed" from the path list by "merging" it and the two branch clusters it joins into a "super cluster". In fact, no actual modification of medial axis voxels occurs, the merging is implemented in a book-keeping manner.

  The cluster merging procedure can in fact be delayed until case 4.12, the fiber tracing procedure, is run.

Output from case4.10 is a trimmed medial axis, which can either be stored in "voxel-location" format files, or preferably in "cluster-path" format which retains all book-keeping information on (merged) clusters and remaining paths. Complete information on the structure of the cluster/path files is found in section 9.4 of the 3DMA General Users Manual.

Test example:
	- download the input file case4.10.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case4.10

	- This will produce the output file case4.10.out and
	  trimmed medial axis data files cp_sngl, cp_loc and 
	  cp_struct of cluster/path format in the ma_t/ directory.

Analysis continues using trimmed medial axis files stored in the directory /nfs/3dma_fiber/Test/ma_t/

The medial axis can be viewed either slicewise (rasterfile format), or in 3D (Geomview or Inventor format). Each format has advantages and disadvantages. The continuity of the medial axis is difficult to follow in slice-by-slice viewing, but the relationship of the medial axis to the grain phase is clearer. 3D viewing is preferable for seeing the continuity of the medial axis, however the simultaneous plotting of the pore-grain surface hides much of the medial axis and makes it somewhat difficult to view the medial-axis/grain surface relationship. Plotting the medial axis in 3D without the grains has the downside that it makes perspective viewing of the medial axis difficult.

We recommend plotting the medial axis in color, where the color spectrum provides information on the burn number distribution along the medial axis. Typically we use a "rainbow" color scale, small burn numbers being red, large burn number being blue/violet.

Test example:
	- download the input file case4.1.in
	  into the /nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case4.1

The case4.1.in file is set up to produce a Geomview plot, ma.list, in the plots/ directory. A view of this 3D plot is shown in Fig. 12.

Fig. 12:  A Geomview graphics view of the 3D medial axis of the sample data set. Grains are not plotted. Note how difficult it is to see perspective in the 3D medial axis plots. Color reflects fiber radius, with green fibers being thicker than red. The radius is digitized in terms of the 7 micron resolution of these images.

Using the relation

	distance = (2*bnum - 0.5)*vox_len

Test example:
	- download the input file case4.4.in
	  into the /nfs/3dma_fiber/Test/cases/ directory.

	% runcase linux case4.4

	- Produces a jgraph plot, results/thick.jgr, of the distribution of (3D)
	  diameters (Fig. 14) for the 10 slice data set.

When the medial axis is computed for a single (2D) slice, the distance measured is the 2D width of the fibers within the plane of the slice. Thus comparing the distribution of these distances for medial axes computed for a representative (2D) slice in the image against that for the entire 3D image can give some indication of the anisotropy in fiber thicknesses. This is especially true for collapsed (ribbon-like) fibers that may lie preferentially in one of the three principal planes of the image.

In order to measure in 2D, a single slice must be extracted from the 3D data stack and a medial axis prepared for the slice. The following example takes slice 5 from the 3D example data stack, and prepares a medial axis for this slice.

Test example:
	- download the input files:
		case2.5slc.in
		case4.1slc.in
		case4.10slc.in
		case4.1slc_t.in
		case4.4slc.in
	  into the /nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case2.5slc

	- Generates the medial axis for slice 5 from the cleaned
	  segmented image in the c_seg directory. The medial axis information
	  for this slice is stored in the ma/fbr_slc.005.gz files and the
	  corresponding burn information in burn/fbr_slc.005.gz

	% runcase linux case4.1slc

	- Produces a raster plot ma/fbr_slc.005.ras.gz (Fig. 13 (left)) of the
	  medial axis for this slice

	% runcase linux case4.10slc

	- Provides some clean-up of the slice medial axis so that 2D widths
	  computed are not overly reflective of small fiber segments in the
	  slice which are really pieces of fibers passing through the slice.
	  The cleaned up medial axis information is stored in the
	  ma_t/slc_{sngl,loc,struct} files.

	% runcase linux case4.1slc_t

	- Produces a raster plot ma_t/fbr_slc.005.ras.gz (Fig. 13 (right)) of
	  the trimmed medial axis for this slice.

	% runcase linux case4.4slc

	- Produces a jgraph plot, results/width.jgr, of the distribution of (2D)
	  widths (Fig. 14) for this slice.

Fig. 13 compares the medial axis obtained for the single slice (slice 5) before and after trimming the medial axis to ignore small fiber segments passing through the slice. Fig. 14 compares the width (one 2D slice) and thickness (a 3D set of 10 slices) distributions obtained. The results tend to support an observation that these fibers have a slight anisotropy in cross sectional diameter. Obviously the test data set is too small, and a larger set of slices, combined with more samples of single slices would be required to confirm this observation.

Fig. 13:  The medial axis obtained for the single slice (slice 5) before (left) and after (right) trimming designed to ignore small fiber segments passing through the slice.

Fig. 14:  Comparison of 2D width and 3D (smallest) diameter of the fibers in the example data set. The values "< X.Y >" denote the measured average diameter of each distribution.

Case 4.12 traces fibers through the 3D medial axis using a minimum bending angle criterion to connect fiber segments across bonding points.

Test example:
	- download the input file case4.12.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% cd /nfs/3dma_fiber/Test/cases/
	% runcase linux case4.12

After fiber reconstruction is finished, case4.12 generates a number of statistics (distribution) on the fibers. All plots of these distributions are now plotted in a single jgraph file (../results/fbr.jgr) on multiple pages. The order of plots are:

	  fiber segment length		(Fig. 15)
	  fiber length			(Fig. 16)
	  fiber curl (tortuosity)	(Fig. 17)
	  bonds per fiber		(Fig. 18)
	  bonds per unit fiber length	(Fig. 18)
	  fiber crossing angle		(Fig. 18)
	  fiber moment of inertia	(Fig. 19)

Fig. 15:  Length distributions (expressed as probability density) of fiber segments for: (top left) all segments; (top right) branch-branch segments; (bottom left) branch-leaf segments; (bottom right) leaf-leaf segments. In this example, all leaf-leaf segments have been removed from the data set by the medial axis trimming case 4.10. The values "(av X.Y)" denote the measured average length of each distribution.

Fig. 16:  Length distributions (expressed as probability density) of fiber (top left) all fibers; (top right) branch-branch fibers; (bottom left) branch-leaf fibers; (bottom right) leaf-leaf fibers. The values "(av X.Y)" denote the measured average length of each distribution.

Fig. 17:  Fiber curl (tortuosity) distributions (expressed as incremental probability): (top left) all fibers; (top right) branch-branch fibers; (bottom left) branch-leaf fibers; (bottom right) leaf-leaf fibers. The values "(av X.Y)" denote the measured average curl of each distribution. The solid line shows the cumulative probability distribution.

Fig. 18:  Distributions for: (left) bonds per fiber; (middle) bonds per unit fiber length; (right) fiber crossing angle. All distributions are expresssed as incremental probability. The bonds per fiber distributions are highly distorted by the small data set size of this example.

Fig. 19:  (Top row) Fiber moment of inertia scatterplots: (top left) smallest MoI vs middle; (top middle) smallest MoI ve largest; (top right) middle MoI ve largest. (Rows 2-4) Spatial distribution of principle axes associated with the: (row 2) smallest MoI; (row 3) middle MoI; (row 4) largest MoI.

        Case 0.1 generates a simulated data set of fibers in the shape of a 3D
        angular fan. To reduce the size of the data set, it exploits spherical
        symmetry (and symmetries of the medial axis algorithm employed by
	3DMA-Fiber).  The fan consists of non-overlapping fibers of size

                  thickness x thickness x length  (units = microns)

        As the zenith angle (phi - measured from the x-y plane) increases,
        fewer azimuthal angles, (theta) are sampled in order to keep the
        fibers from overlapping.
 
        The zenith angle, phi, is broken into D_PHI degree increments.
	(D_PHI is a defined value currently set to 1 degree.)
 
        For a given phi, the angle between the vectors
            v1 = (cos(phi),0,sin(phi))
        and
            v2 = (cos(phi),cos(phi),sqrt(2)sin(phi))/sqrt(2) (theta = 45 degrees)
        is
	    alpha(phi) = arccos(v1 dot v2).

        The angle alpha(phi) is broken into D_ALPHA degree increments.
	(D_ALPHA is a defined value currently set to 1 degree.)
 
        The array of fibers generated is not "rectangular" in phi,alpha.
	Taking into account the symmetries of the circle, the range of
	(theta,phi) angles swept out by the fiber array is:
 
        for fixed theta in 0 <= theta <= pi/4
        then               0 <=   phi <= atan(cos(theta)) ]
 
        Written for fixed phi:
 
        for 0 <= phi <= atan(1/sqrt(2))    then 0 <= theta <= pi/4
        for atan(1/sqrt(2)) <= phi <= pi/4 then 0 <= theta <= acos(tan(phi))

A simulated data set consists of straight line segments, corresponding to the start and end points of the center line of "free fiber segments". A free fiber segment is a segment of a fiber between bonding (i.e. touching) points with other fibers, or between the start (end) of a fiber and its first bonding point. The simulated fiber data set consists of a simple ascii printout of the free segments, with a header consisting of global fiber information. The precise format of a simulated data file is

	float - fiber thickness (microns assumed)
	float - fiber width (microns assumed)
	float - fiber length (microns assumed)
	float - x-offset to origin
	float - y-offset to origin
	float - z-offset to origin
	int   - number of free fiber segments to follow
	6 floats - x y z coordinates of start point, x y z coordinates of end point of segment 1.
	6 floats - x y z coordinates of start point, x y z coordinates of end point of segment 2.
	  ...

Notes:

• The fibers are generated in a coordinate system (whose units are assumed to be microns). All fibers lie in a rectangular cube (xmin,xmax) x (ymin,ymax) x (zmin,zmax). The three offset values stored in the header are respectively xmin, ymin and zmin.

Our discussion of case0.1 will not only consider generation of the fan simulated data set, but will illustrate the entire processing sequence through to the generation of the length correction table.

  cd /nfs/3dma_fiber
  mkdir fan

- Organize the following subdirectories under /nfs/3dma_fiber/fan/:

   cd Test
   mkdir burn_fan cases_fan ma_fan ma_fan_t plots_fan results_fan seg_fan tom_fan

- Download the archived set of fanX.Y.in input files into the /nfs//3dma_fiber/fan/cases_fan subdirectory.

- Extract the input files.
  unzip fbr_fan.tar.gz
  tar -xvf fbr_fan.tar

Fan generation and length correction table example:

	% cd /nfs/3dma_fiber/fan/cases_fan/
	% runcase linux fan0.1
	  - Generates a fan simulated data set "f_fan" in the cases_fan
	    directory
	% runcase linux fan0.3
	  - Produces a 3D digitized greyscale (tomographic) data representation
	    of the fan at 4 micron voxel size ("resolution"). The tomographic
	    data is stored in the tom_fan directory.
	% runcase linux fan1.1
	  - Produces rasterfile plots of the data set in the tom-fan directory
	    for viewing with "xv".
	% runcase linux fan1.4
	  - Produces a histogram of greyscale intensities in the digitized data
	    set for viewing via "jgraph". Provides guidance for choose
	    thresholding paramenters in kriging segmentation step.
	% runcase linux fan1.7
	  - Segmentation step via indicator krigint. Segmented data stored in
	    seg_fan directory.
	% runcase linux fan2.3
	  - Analysis of disconnected volumes in the segmented image. For the
	    fan data, the volume of each disconnected "pore" actually
	    corresponds to a volume (in voxels) of a digitized fiber in the fan.
	% runcase linux fan2.5
	  - Generates medial axis and burn data for each fiber.
	    See Notes below.
	    Stores burn data in burn_fan directory and medial axis data in
	    ma_fan directory.
	% runcase linux fan4.9
	  - Converts medial axis data from voxel list format to cluster/path
	    format.
	% runcase linux fan4.1
	  - Generates a 3D view of the medial axis of the fiber fan for
	    viewing by geomview. (See Fig. 20.) The plot is put in the plots_fan
	    directory.
	% runcase linux fan4.11
	  - Generate path length and coordination number distributions of
	    medial axis. Provides information for medial axis trimming step.
	% runcase linux fan4.10
	  - Trims small branch-leaf paths from the fan medial axis. These
	    small protrubances result from the medial axis algorithm respecting
	    roughness on the surface of the digitized fibers. Trimmed medial
	    axis data stored in the ma_t_fan directory.
	% runcase linux fan4.12
	  - Identify fibers in the fan from the corrected medial axis and
	    outputs a length correction table, ratio_table, in the cases_fan
	    directory.
	    See Notes below.
	% runcase linux fan_corr4.12
	  - Applies the length correction table to the fan as a check.
	    See Notes below.

Fig. 20:  Geomview view of the medial axis of the simulated fiber fan.

Notes:

• The medial axis algorithm used in the 3DMA-Fiber code (Lee, Kashyap and Chu. See in References in 3dma_fiber.html) has known orientation biases. While, intuititely, the medial axis algorithm reduces the thickness of a fiber (down to a central axis) it can also reduce the length of a fiber. The amount of reduction depends on the curvature of the fiber, its tapering near the ends, and its angular orientation in space. The worst case is for straight, digitized, constant thickness fibers oriented at the angle theta = 0, phi = 45 degrees. At this angle the legnth reduction can be substantial.
• The length correction table consists of a multiplicative correction factor computed from the fan for each (phi,alpha) angular increment (currently set to 1 degree in each angle). The output table is in binary, however its contents, in ASCII, are echoed into the case4.12.out file. Statistics on the length correction factors are also printed in the case4.12.out. In addition to computing the mean and standard deviation value of the correction factors, a quartile based "box-whisker" analysis is done, giving values for the q_25, q_50 and q_75 correction factors, as well as the number of correction factors in the table that at mild of extreme outliers.
• One would therefore expect that application of the length correction table to the digitized fan data (which after-all is used to create the correction table) would lead to exact legnth corrections for all fibers in the fan. However, a subtlety occurs. Creating the length correction table (1 entry per fber in the fan) requires identifying the angular orientation of each digitized fiber medial axis. This angular orientation is determined by the angle of the vector between the origin and the end point (closest to the origin) of the fiber's medial axis. This method has been determined to be robust. However, when the length correction table is applied to correct a general fiber (i.e. associated with a users data set and not any specific fan) the angle of the fiber should be determined by its end-to-end orientation in space (computed by moment of inertia analysis using the principle direction associated with the smallest moment). FOR THE MEDIAL AXIS OF FAN FIBERS, the origin-to-end-point angle and the moment of inertial principle direction may in fact slightly disagree (details depend on the medial axis construction for each particular oriented fiber). Thus application of the length correction table to the fiber fan will not, in general, produce an exact correction.

Case 0.2 produces a "trellised" fiber fan. It is based upon the fiber fan generated by case0.1, but keeps only fibers that are 5 degrees apart (in phi, alpha). In addition, it runs two fibers perpendicularly across those fibers having the same phi-value. Thus fibers at each phi-value are bound together in a small network. This data set is provided to test the fiber finding algorithms of 3DMA-Fiber. Our discussion of case0.2 considers generation of the trellis simulated data set, as well as the entire processing sequence through to fiber identificaiton.

  cd /nfs/3dma_fiber
  mkdir trellis

- Organize the following subdirectories under /nfs/3dma_fiber/fan/:

   cd Test
   mkdir burn_trellis cases_trellis ma_trellis plots_trellis results_trellis seg_trellis tom_trellis

- Download the archived set of treX.Y.in input files into the /nfs//3dma_fiber/fan/cases_trellis subdirectory.

- Extract the input files.
  unzip fbr_trellis.tar.gz
  tar -xvf fbr_trellis.tar

Trellis generation and analysis example:

	% cd /nfs/3dma_fiber/trellis/cases_trellis/
	% runcase linux tre0.2
	  - Generates a trellis simulated data set "f_tre" in the cases_trellis
	    directory
	% runcase linux tre0.3
	  - Produces a 3D digitized greyscale (tomographic) data representation
	    of the trellis at 4 micron voxel size ("resolution"). The
	    tomographic data is stored in the tom_trellis directory.
	% runcase linux tre1.1
	  - Produces rasterfile plots of the data set in the tom_trellis
	    directory for viewing with "xv".
	% runcase linux tre1.4
	  - Produces a histogram of greyscale intensities in the digitized data
	    set for viewing via "jgraph". Provides guidance for choose
	    thresholding paramenters in kriging segmentation step.
	% runcase linux tre1.7
	  - Segmentation step via indicator krigint. Segmented data stored in
	    seg_trellis directory.
	% runcase linux tre2.3
	  - Analysis of disconnected volumes in the segmented image. For the
	    trellis data, the volume of each disconnected "pore" actually
	    corresponds to a volume (in voxels) of the digitized fiber trellis
	    at each angle phi.
	% runcase linux tre2.5
	  - Generates medial axis and burn data for each fiber.
	    Stores burn data in burn_trellis directory and medial axis data in
	    ma_trellis directory.
	% runcase linux tre4.9
	  - Converts medial axis data from voxel list format to cluster/path
	    format.
	% runcase linux tre4.1
	  - Generates a 3D view of the medial axis of the fiber fan for
	    viewing by geomview. (See Fig. 21.) The plot is put in the
	    plots_trellis directory.
	% runcase linux tre4.11
	  - Generate path length and coordination number distributions of
	    medial axis.
	% runcase linux tre4.12
	  - Identify fibers in the trellis from the medial axis and
	    apply the length correction table, ratio_table, in the cases_fan
	    directory.
	    
	    YOU WILL HAVE TO EDIT THE tre4.12.in FILE AND SET THE PATH TO
	    THE LENGTH CORRECTION TABLE APPROPRIATE FOR YOUR DIRECTORY TREE.
	    SEE THE LINE

	    Enter basename for input ratio table file(s):

	    in tre4.12.in.  IF YOU ARE FOLLOWING THE EXMPLES HERE, THE
	    APPROPRIATE ANSWER SHOULD BE

	    Enter basename for input ratio table file(s): ../../fan/cases_fan/ratio_table
	    

Fig. 21:  Geomview view of the medial axis of the simulated fiber trellis.

Case 0.3 converts a simulated fiber data set to a digitized greyscale image. Recall (see notes under case0.1) that a simulated fiber set consists of a list of (the end point coordinates of) straight center-lines of free fiber segments. Case0.3 digitizes each free fiber segment according to the thickness and width specified for the fibers in the header information and according to the voxel size specified by the user.

For example, a fiber segment that is 100 microns long, and represents a fiber that is 20 microns in thickness and 30 microns in width, would be digitized at 4 micron resolution, as a rectangle occupying 5 x 6 x 20 voxels.

This simple example assumes the fiber is to align to the coordinate directions of the digitized space. If is lies at angles to the coordinate directions, the digitiztion of the free segment must take this into account in assigning "partial voxels". The assignemnt of partial voxels is done using greyscale values. Weighting the fraction of a voxel occupied by an appriate fraction of maximum grey-scale value.

The details of the digitization are the following.

• Each free fiber segment center line is rotated/translated to lie along the x-axis, with one end at the origin.
• The free fiber segment is digitized into voxels aligning with the coordinate axis.
• Each voxel of the aligned, digitized fiber segment is subdived into n x n x n subvoxels ("n" is specified by the user by specifying the linear dimension of a subvoxel)
• Each subvoxel is accorded n^(-3) of the maximum grey-scale value.
• The center of each subvoxel is then rotated/translated back to the appropriate direction/position in 3-space. The voxel in which this transformed center-point lies accumulates an n^(-3) amount of greyscale value to indicate partial volume occupation by this piece of the digitized fiber segment.

No specific input example is provided here. Examples of case0.3 runs are given in the discussion of case0.1 and case0.2.

3DMA-Fiber has the ability to compute the tortuosity through the void space. To illustrate this capability, we provide a different test sample consisting of a 70 x 70 x 70 voxel cubic subvolume of the unwoven fiber mesh. (For comparing tortuosities in different directions, it is preferable that sample lengths be the same in all directions.)

Directory set-up, download example input files

For clear organization, we create new subdirectories in the /nfs/3dma_fiber/Test/ directory for the void space analysis.

   cd Test
   mkdir v_burn v_cases v_c_seg v_ma v_ma_t v_plots v_results

These directories will store files analogous to those in the preceeding fiber example.

For this example we provide a data set that has already been, resized (case1.2), segmented (case1.6), and had disconnected volumes identified and removed (cases 2.3, 2.4).

Use the following link to download the archived (fbr1-cseg.tar) data set of example segmented input files into the /nfs/3dma_fiber/Test/v_c_seg/ directory.

This archived data set is 82 KBytes in size.

The data set contains one segmented volume file, 'fbr1.gz', containing a 70 x 70 x 70 volume of the same unwoven fiber mat sample analyzed above. The data are stored in the bit-packed-byte format we use for segmented files. The format for this file is given in section 9.1 of the 3DMA General Users Manual. Except for 4 leading integers, the rest of the data is bytes (unsigned chars) and is therefore independent of machine endian. The leading 4 integers are, of course, endian dependent and are stored in little-endian format. If you are using a big-endian machine, do the appropriate byte swapping of the first 4 integers in the data file.

- Extract the segmented files.
  unzip fbr1-cseg.tar.gz
  tar -xvf fbr1-cseg.tar

Process the data in preparation for the tortuosity computation

Test example:
	- download the input files:
		case2.5v.in
		case4.10v.in
		case4.1vt.in
		case4.9v.in
	  into the /nfs/3dma_fiber/Test/v_cases/ directory.

	% cd /nfs/3dma_fiber/Test/v_cases/
	% runcase linux case2.5v

	- creates the medial axis for the VOID space of the fiber mesh.
	  (Note in particular the

		Invert data (y,n(dflt))?: n

	  line in case2.5v.in which guarantees the medial axis is computed
	  for the void space.) Puts the medial axis information in the
	  v_ma/fbr1.xxx.gz files.

	% runcase linux case4.10v

	- trims the medial axis. Some hint on length parameters used to control
	  medial axis trimming can be gained by running Case 4.11 (histogram
	  MA info) and requesting the histogram of path lengths. Download the
	  input file case4.11v.in for an example.  Puts the trimmed medial axis
	  information in the v_ma_t/cp_{sngl,log,struct) files.

	% runcase linux case4.1vt

	- Prepares a geomview plot (../v_plots/v_ma_t.list) of the trimmed,
	  void space medial axis (Fig. 22)

	% runcase linux case4.9v

	- Converts the cluster/parth formated v_ma_t/cp_{sngl,log,struct) files
	  to voxel location list format (v_ma_t/fbr1.xxx.gz) since the case 4.6
	  runs to be performed next have not been upgraded (yet) to input the
	  cluster path format.

Fig. 22:  A Geomview plot of the trimmed void space medial axis for the 70x70x70 voxel data set.

Tortuosity computation

With a trimmed medial axis computed for the pore space inside a rectangular volume and output in MA voxel list format, tortuosity distributions for paths between any two faces of the rectangular volume can be computed. Here is a quick description of the algorithm used to compute geometrical tortuosities.

Fig. 23:  Labeling of 3D image faces in Case 4.6.

For each MA voxel i on face 3 and for each MA voxel j on face 4, the algorithm will find the shortest path (if one exists) from i to j through the medial axis network. Let l_ij be the length of this shortest path. and d_34 be the Euclidean distance between faces 3 and 4. The geometrical tortuosity t_ij of the path is defined as t_ij = l_ij/d_34. Geometrical tortuosities are computed for the shortest paths between all possible entrance/exit voxel pairs i, j on the MA network.

Tortuosity results are summarized in two ways. The first is the geometrical tortuosity distribution for all entrance/exit voxel pairs:

   {t_ij   |  i is an entrance (face 3) MA voxel, j is an exit (face 4) MA voxel},

The second is the geometrical tortuosity distribution for the shortest path connecting each entrance voxel i to the exit face:

   {t_i = min_j (t_ij)   |  j is an exit face MA voxel},

where i ranges over all MA voxels on the extrance face.

Test example:
	- download the input file case4.6vx.in into the
	/nfs/3dma_fiber/Test/cases/ directory.

	% runcase linux case4.6vx

	- This outputs the Geomview plots vtort_x_all and vtort_x_abs in
	  ../v_plots, and the distribution graph vtort_x.jgr in ../v_results.

	% cd ../v_plots
	% geomview tort_x_all 

	- This will plot a set of paths, each path being the shortest path
	  joining each i,j (entrance face 3, exit face 4) voxel pair. This plot
	  is shown top left in Fig. 24.

	% geomview tort_x_abs 

	- This will plot a set of paths, each path being the shortest path
	  between each i (entrance face 3) voxel and the exit face 4. This plot
	  is shown top right in Fig. 24.

	% cd ../v_results 
	% jgraph < vtort_x.jgr > vtort_x.ps
	% ggv vtort_x.ps

	- This plots the tortuosity distributions for the paths in tort_x_all
	  (top plot) and tort_x_abs (bottom plot). These two plots are shown
	  on the left in Fig. 25.

	- download the input files
		case4.6vy.in
		case4.6vz.in
	into the /nfs/3dma_fiber/Test/cases/ directory.

	% runcase linux case4.6vy
	% runcase linux case4.6vz

	- These output analogous data for path tortuosities connecting the
	  y- and z-direction faces. These paths and their tortuosity
	  distributions are also shown in Figs. 24 and 25.

Fig. 24:  Geomview plots of the medial axis paths from: (top) faces 3 -> 4 (x-direction); (middle) faces 5 -> 6 (y-direction); (bottom) faces 1 -> 2 (z-direction). For each direction, the left-hand plot shows the shortest paths between every inlet-outlet voxel pair while the right-hand plot shows only the shortest paths across the image from each inlet voxel. Note that the y-direction points toward the viewer.

Fig. 25:  Distributions of void space tortuosities in the (left) x- (middle) y- and (right) z-directions from the plots of Fig. 24. Top histograms correspond to the plots on the left of Fig. 24, bottom histograms to the plots on the right of Fig. 24.

In an analogous manner, starting with the segmented fiber volume c_seg/fbr.gz, one can also compute paths through the medial axis of the fiber nextwork and obtain fiber network tortuosities. (The changes to the files case2.5v.in through case4.6zv.in needed to do this should be straightforward.)
Note: this is NOT the same tortuosity as the fiber curl statistics computed in case4.12.
In Fig. 26, we show the tortuosity distributions we obtain for the fiber network for this data set. Note the higher tortuosity of fiber paths in the z-direction. In contrast to the void space, the x-y layering of the fibers makes it difficult to find "straight" paths through the fiber mesh in the z-direction.

Fig. 26:  Distributions of fiber network tortuosities in the (left) x- (middle) y- and (right) z-directions from the c_seg/fbr.gz data files.

• 3DMA General Users Manual,  W.B. Lindquist.
• Three-dimensional image analysis of fibrous material,  H. Yang and W.B. Lindquist
• A geometric analysis of 3D fiber networks from high resolution images,  H. Yang and W.B. Lindquist