// 20150902 Akira Muto
// Motivation: This macro was written to automate counting of the paramecia left in the recording chamber in a time-lapse movie, in which a zebrafish larva eats the paramecia over time.
// Each experiment results in an AVI file (1024x1024 pixesl x 6600 frames, 6.4GB). Dozens of recordings in a day. 
// Processing on a movie file takes about 20 min. Processing of 1-day recordings takes one entire day or even more. Better to automate it. 
//
// What it does: My desktop PC has 32GB memory but it is not enough to process (image division) one file (6600 frames) and I need split each file into 2 (3300 frames and 3300 frames) to avoid the memory error.
// For background subtractio, the movie is divided by the averaged image pixel-by-pixel. The resulted ratio images are saved as an AVI file.
// The ratio images are thresholded at the pixel value of 94 to separate the paramecia from the background. 
// Because the recording chamber is circular, an Oval ROI is set for the particle analysis. 
// Paramecia are detected by Particle Analysis with the parameter set of size=18-72 circularity=0.5-1. 
// The particle count is saved as a text file in the same folder. 


//////////////// Specify the image folder
source_dir = getDirectory("Please select the data folder.");
print("Data folder:"+source_dir);

//////////////// Obtain the names of the files in the folder. 
imagelist = getFileList(source_dir);


////////////////// Batch processing to work on all the AVI files in the folder     
for(i=0; i<imagelist.length; i++){

	if(isImage(imagelist[i])){ // To work on only AVI files. 

	print("imageList["+i+"]: "+imagelist[i]);
	imagefile=replace(source_dir+imagelist[i],File.separator,File.separator);//Dummy to get a string
	run("AVI...", "select=imagefile  first=1 last=6600");
	
	if(nSlices() > 6600) numOfSlice=6600; // To avoid out of memory error, limit to 6600 frames at most. 
	else numOfSlice=nSlices();
	print("numOfSlice="+numOfSlice);

	run("Z Project...", "start=1 stop=numOfSlice projection=[Average Intensity]");
	selectWindow(imagelist[i]);


	/////////// Separte the movie into halves by duplication to avoid the memory error during image division.

	numOfSliceHalf=round(numOfSlice/2);
	initialSliceOfSecondHalf=numOfSliceHalf+1;

	run("Duplicate...", "duplicate range=1-numOfSliceHalf");
	firstHalf = getImageID();
	selectWindow(imagelist[i]);
	
	run("Duplicate...", "duplicate range=initialSliceOfSecondHalf - numOfSlice");
	secondHalf = getImageID();

	selectWindow(imagelist[i]);close();


	///////////////////// Process on each half. Images are divided by the averaged image for background subtraction. 

	imageCalculator("Divide create 32-bit stack", firstHalf,"AVG_"+imagelist[i]);	
	run("Enhance Contrast", "saturated=0.35");
	run("8-bit");
	selectImage(firstHalf);
	close();
	
	imageCalculator("Divide create 32-bit stack", secondHalf,"AVG_"+imagelist[i]);	
	run("Enhance Contrast", "saturated=0.35");
	run("8-bit");
	selectImage(secondHalf);
	close();

	selectWindow("AVG_"+imagelist[i]);
	close();

	
	//////////////// Concatenate the two halves. 

	firstResult="Result of "+firstHalf;
	secondResult="Result of "+secondHalf;
	//ConcatenatedStack=replace(imagelist[i],".avi","_ratio.avi");//Somehow large-sized (>6GB) AVI file was not saved properly.
	ConcatenatedStack=replace(imagelist[i],".avi","_ratio.raw");// Even a >6GB stack was successfully saved (only) in Raw format.

	run("Concatenate...", "  title=[ConcatenatedStack] image1=[firstResult] image2=[secondResult] image3=[-- None --]");


	/////////// Save the resulted ratio images.

	pathtosave=source_dir+ConcatenatedStack;
	//run("AVI... ", "compression=Uncompressed frame=30 save=pathtosave");//Somehow large-sized (>6GB) AVI file was not saved properly.
	//saveAs("Tiff", pathtosave);//Strangely, large-sized TIFF (>6GB) could not be saved (at all). 
	saveAs("Raw data",pathtosave);// Even a >6GB stack was successfully saved (only) in Raw format.
	print("Ratio image was saved as "+pathtosave);


	////////// Particle analysis to detect paramecia 

	setThreshold(0, 140);//Determine for your own movie.e.g., 94 or 140. 
	run("Convert to Mask", "method=Default background=Light");
	run("Invert", "stack");
	makeOval(0, 0, 1023, 1023);
	run("Analyze Particles...", "size=18-72 circularity=0.5-1 show=Masks display exclude clear summarize stack");


	////// Save the result of the particle analysis

	SummaryWindowName="Summary of "+ConcatenatedStack;
	imagefilename=source_dir+imagelist[i];
	//summaryfilename=replace(source_dir,".avi","_summary.txt");
	summaryfilename=substring(imagefilename,0,lengthOf(imagefilename)-4);//to remove ".avi". avi appears in behavior (folder name!)
	selectWindow(SummaryWindowName); 
	saveAs("Text", summaryfilename+"_summary.text");
	run("Close");

	print("Summary was saved as "+summaryfilename+"_summary.text");

	
	///////////// Close the image windows and the Results window. 

	while (nImages>0) {selectImage(nImages); close();}//Closing all the image windows 
	selectWindow("Results");run("Close");
	 

	}// End of if-clause

}// End of for-loop




////////////// This function was borrowed from http://adm.mri.cnrs.fr/datas/fichiers/articles/60/184.pdf 
function isImage(filename) {
 extensions = newArray("avi");
 result = false;
 for (i=0; i<extensions.length; i++) {
 if (endsWith(toLowerCase(filename), "." + extensions[i])) 
result = true;
 }
 return result;
}