Biofilm thickness (z-stacks)

I have 3D images from confocal imaging of my biofilms. they are taken through z stacks. now I want to measure the thickness of biofilm I have. How I can do that?

Hi Sababiba,

You can copy and paste the macro code below in the script editor (and select Language IJ1 macro ).

//open an image 
run("Confocal Series (2.2MB)");
image_Name = getTitle();

//get some informations about the image
getDimensions(image_width, image_height, image_channels, image_slices, image_frames);
getVoxelSize(voxel_width, voxel_height, voxel_depth, voxel_unit);

// here we'll make a line in the middle of the image
y_line = image_height/2;
makeLine(0, y_line ,image_width , y_line);

// then we'll create a new iamge by re-slicing
// re-sampling through this line in the z-direction
run("Reslice [/]...", "output="+voxel_depth+" slice_count=1");
selectWindow("Reslice of "+image_Name);

// Now you could draw a line and get the length of this line 
// ther is a live measure in the main ImageJ/Fiji bar when you draw it 
// or you can measure
setTool("line");
while (roiManager("Count") < 1  ){
	waitForUser("Please make a line");
		if (selectionType > -1 ){
		roiManager("Add");
		roiManager("Measure");
	}
}

BUT this manual way even if it’s quick can lead to some bias in your measurements, because You have to decide WHERE you place the begin and the end of the line ?

When you want to measure the diameter of point object like a bead, what is commonly accepted is to measure the Full width at half maximum (have look at FWHM on wikipedia) using a Gaussian fit, but you may need to use a Super Gaussian fit , see below.

// do some cleaning of the environment before starting
run("Close All");
run("Clear Results");

// open an image 
// we need a sampling that is larger than the object
// in the previous image stop there were not enough z 
// above and below the object
run("T1 Head (2.4M, 16-bits)");
image_Name = getTitle();

//get some informations about the image
getDimensions(image_width, image_height, image_channels, image_slices, image_frames);
getVoxelSize(voxel_width, voxel_height, voxel_depth, voxel_unit);

// here we make a line
makeLine(0, 90, 256, 90);

// then we create a new image by re-slicing
run("Reslice [/]...", "output="+voxel_depth+" slice_count=1");
selectImage(nImages);
resliceImage = getTitle();

//Let's make a line
makeLine(25, 67, 240, 67);
waitForUser("We'll make a profile through this line and 'fit' a gaussian formula on it");
fitGaussian();
waitForUser("The inhomogenity of the signal causes trouble ");

// Let's blur a bit the image
selectImage(resliceImage);
run("Duplicate...", "title=["+resliceImage+"-Blurred]");
run("Gaussian Blur...", "sigma=5");
waitForUser("So we blur the image a bit before doing the fitting");
makeLine(25, 67, 240, 67);
fitGaussian();

selectImage(resliceImage+"-Blurred");
waitForUser("But this is still not 'perfect' \nbecause the length of the line you used \nwill affect your fitting and the results");
makeLine(0, 67, image_width, 67);
fitGaussian();



selectImage(resliceImage+"-Blurred");
waitForUser("One possible way is to fit a 'SuperGaussian', \nbut it needs estimates(see function below)");
makeLine(0, 67, image_width, 67);
fitSUPERGaussian();


// companion function to make the math for FWHM
function fitGaussian(){
	
	y = getProfile();
	x = Array.getSequence(lengthOf(y));
	Fit.doFit("Gaussian", x, y) ;
	Fit.plot;
			
	sortedParameter = Fit.p(3); // parameter d of gaussian
	rSquared = Fit.rSquared ; 
			
	FWHM = (2 * sqrt( 2 * log(2) ) ) * sortedParameter ;// http://fr.wikipedia.org/wiki/Largeur_%C3%A0_mi-hauteur
	setResult("FWHM ("+voxel_unit+")", nResults, FWHM * voxel_height);
	setResult("Label",nResults-1,image_Name);
	setResult("rSquared",nResults-1,rSquared);
	updateResults();
	selectWindow("Results");
}


function fitSUPERGaussian(){

	getVoxelSize(voxel_width, voxel_height, voxel_depth, voxel_unit);

	yCurrent = getProfile();											
	xPlot = Array.getSequence(lengthOf(yCurrent));			
			
	pidiv2 = PI/2;
	formulaString = "y = a + b * exp( "+pidiv2+"*(-pow ( abs(x-c), e/2 ) ) / ( pow( d* "+pidiv2+" , e/2 ) ) )";
	
	pourcentOfMax 				= 50	; // percent of the Maxmimal value of the Gaussian you want to consider as the diameter, 50 = FWHM
	fractionOfMax 				= 1/pourcentOfMax *100; 
	
	baseline=0;
	maxOfCurve = 255;
	centerOfCurve = 125 ;
	fullWithFactor = 200; // estimated diameter
	gamma = 2 // is close to 2 for a gaussian
	initialGuesses = newArray(baseline, maxOfCurve, centerOfCurve, fullWithFactor, gamma);
	
	Fit.doFit(formulaString, xPlot, yCurrent,initialGuesses);
	Fit.plot();	
	parameterA	= Fit.p(0); 	// baseLine 
	parameterB 	= Fit.p(1);	// top - baseLine
	parameterC  = Fit.p(2);		// center
	parameterD 	= Fit.p(3);	// function of FW 
	parameterE 	= Fit.p(4);	// N, the great lord gamma
	rSquared 	= Fit.rSquared;
	
	radius 	=  parameterD * pow( (2/PI), (2/parameterE) - 1 ) * pow ( log(fractionOfMax), 2/parameterE) ; // in pixel
	
	setResult("FWHM ("+voxel_unit+")", nResults, 2*radius*voxel_width);
	setResult("Label",nResults-1,image_Name);
	setResult("rSquared",nResults-1,rSquared);
	updateResults();
	selectWindow("Results");

}

I hope it will help you to go forward.

Cheers,

Romain

1 Like