// @File(label="Input directory", description="Select the directory with input images", style="directory") inPath
// @String(label="Scale Factor", description="Scale to this fraction. Runs faster but less accurate", choices={"1","0.5","0.25"}, default="0.25",persist=false ) scaleFrac
// @String(label="Verbose Level", description="Higher number prints more info", choices={"0","1","2"}, default="0",persist=false ) verbose

//
// --------------------------------------------------
// Analysis of Elastin in Aorta sections
// --------------------------------------------------
//
// Script which analyses pre-segmented images by skeletonisation.
// Requires images pre-processed with Ilastik to produce Simple Segmentation masks
// in 8-bit tif format.
//
// Regions labelled '1' will be used for further processing, BG can be anything else.
//
//												
//									- Dave Mason [dnmason@liv.ac.uk]
//										Centre for Cell Imaging, Liverpool [http://cci.liv.ac.uk]
//
//---------------------------------------------------

//-- RUN OPTIONS:

//------------------------------------
//-- Define verbosity:
//	0= Finished message only
//	1= File Processing (show images, give updates on processing)
//	2= Full debug output (as above but also lists arrays)
verbose=parseInt(verbose);
//verbose=0;
//------------------------------------

if (verbose==0){setBatchMode(true);}

//------------------------------------
//-- When refining the mask, only include objects bigger than this (pixels)
minObjSize=80;
//------------------------------------

//--------------------------------//------------------------------------
//-- PART 1: Get file list and specify arrays
//--------------------------------//------------------------------------
saveSettings();
if (verbose>0){print("------------------------");}
if (verbose>0){print("Processing "+inPath);}
dirListAll = getFileList(inPath);

if (verbose>0){print("Found "+dirListAll.length+" total files");}
if (verbose>1){Array.print(dirListAll);}
//-- Refine dirList to only include requested filetypes
dirList=newArray("removeMe");
for (q = 0; q < dirListAll.length; q++) {
	if(endsWith(dirListAll[q], "tif")){
			dirList=Array.concat(dirListAll[q],dirList);
			if (verbose>1){print("Including "+dirListAll[q]);}
	}else{
			if (verbose>1){print("Excluding "+dirListAll[q]);}
		}
}


//-- Sanity check, dirList will always have at least one entry
if (dirList.length==1){
	//-- No eligible files found
	if (verbose>0){
		print("[ERROR] No eligible files found. Quitting.");
		exit("No eligible files found!")
				}	
	}else{
	//-- Remove the dummy line and continue
	dirList=Array.trim(dirList,dirList.length-1);
}
if (verbose>0){print("Found "+dirList.length+" eligible files");}
if (verbose>1){Array.print(dirList);}
//breakplease

//-- Register arrays here with the length equal to dirList
arrCount=newArray(dirList.length);
arrNumBranchAV=newArray(dirList.length);
arrNumBranchSD=newArray(dirList.length);
arrBranchLenAV=newArray(dirList.length);
arrBranchLenSD=newArray(dirList.length);
arrMaxBranchLenAV=newArray(dirList.length);
arrMaxBranchLenSD=newArray(dirList.length);
arrLSP_AV=newArray(dirList.length);
arrLSP_SD=newArray(dirList.length);
//-- and a second set for the cropped regions
arrCount_crop=newArray(dirList.length);
arrNumBranchAV_crop=newArray(dirList.length);
arrNumBranchSD_crop=newArray(dirList.length);
arrBranchLenAV_crop=newArray(dirList.length);
arrBranchLenSD_crop=newArray(dirList.length);
arrMaxBranchLenAV_crop=newArray(dirList.length);
arrMaxBranchLenSD_crop=newArray(dirList.length);
arrLSP_AV_crop=newArray(dirList.length);
arrLSP_SD_crop=newArray(dirList.length);


//--------------------------------//------------------------------------
//-- PART 2: Loop through file list
//--------------------------------//------------------------------------
for (i=0; i<dirList.length; i++){

//-- check requirements for brackets here (to deal with filenames containing spaces)
open(inPath+File.separator+dirList[i]);

if (verbose>0){
print("------------------------");
print("Processing file ["+(i+1)+"/"+dirList.length+"]: "+inPath+File.separator+dirList[i]);}

title=getTitle();

//-- Check whether the image is horizontal or vertical and rotate as necessary
getDimensions(w, h, c, s, f);
if (w>h){
	if (verbose>0){print("Image is landscape, rotating");}
	run("Rotate 90 Degrees Right");
	getDimensions(w, h, c, s, f);
}
if (w<h){
	if (verbose>0){print("Image is portrait");}
}

//-- Convert from a count mask (Ilastik Output) to a binary mask
setThreshold(1, 1);
setOption("BlackBackground", true);
run("Convert to Mask");

//-- Refine mask on full size image
run("Erode");
run("Analyze Particles...", "size="+minObjSize+"-Infinity pixel show=Masks in_situ");

//-- [Optionally] Scale the image to speed up runtime
if (scaleFrac<1){
run("Scale...", "x="+scaleFrac+" y="+scaleFrac+" interpolation=None create title=scaled");
close(title);
selectWindow("scaled");
rename(title);
//-- Refresh with new dimensions so the medial crop works
getDimensions(w, h, c, s, f);
}

//-- Duplicate the mask to run total vs central analysis
run("Duplicate...", "title=medial");
medialFrac=0.15;
boxHeight=floor(h*medialFrac);
run("Specify...", "width="+w+" height="+boxHeight+" x=0 y="+(floor(h/2)-floor(boxHeight/2) ) );
run("Crop");

//-- set up a loop to run the analysis on the whole image (j==0) or the medial subset (j==1)
for (j=0;j<2;j++){

if (j==0){selectWindow(title);if (verbose>0){print("Processing whole image...");}}
if (j==1){selectWindow("medial");if (verbose>0){print("Processing cropped image...");}}

run("Skeletonize");
run("Analyze Skeleton (2D/3D)", "prune=none calculate");

//-- Add the summary stats to the end of the results table
run("Summarize");

//-- catch zero results in first image (no results table) and no skeletons found
//if (isOpen("Results") | nResults>0) {

//-- Record the stats into an array which can be saved at the end as a summary
if (j==0){

	//-- Check that there are more than 1 segment (with only one or zero, summarise doesn't work)
	if (nResults>1){
arrCount[i]=nResults-4;
arrNumBranchAV[i]=getResult("# Branches",nResults-4); //-- Mean
arrNumBranchSD[i]=getResult("# Branches",nResults-3); //-- SD
arrBranchLenAV[i]=getResult("Average Branch Length",nResults-4)/scaleFrac; //-- Mean
arrBranchLenSD[i]=getResult("Average Branch Length",nResults-3)/scaleFrac; //-- SD
arrMaxBranchLenAV[i]=getResult("Maximum Branch Length",nResults-4)/scaleFrac; //-- Mean
arrMaxBranchLenSD[i]=getResult("Maximum Branch Length",nResults-3)/scaleFrac; //-- SD
arrLSP_AV[i]=getResult("Longest Shortest Path",nResults-4)/scaleFrac; //-- Mean
arrLSP_SD[i]=getResult("Longest Shortest Path",nResults-3)/scaleFrac; //-- SD
if (verbose>0){print("Number of Results: "+arrCount[i]);}
	} //-- nResults>0

}
if (j==1){

	//-- Check for no-segments here (empty results table)
	if (nResults>0){
arrCount_crop[i]=nResults-4;
arrNumBranchAV_crop[i]=getResult("# Branches",nResults-4); //-- Mean
arrNumBranchSD_crop[i]=getResult("# Branches",nResults-3); //-- SD
arrBranchLenAV_crop[i]=getResult("Average Branch Length",nResults-4)/scaleFrac; //-- Mean
arrBranchLenSD_crop[i]=getResult("Average Branch Length",nResults-3)/scaleFrac; //-- SD
arrMaxBranchLenAV_crop[i]=getResult("Maximum Branch Length",nResults-4)/scaleFrac; //-- Mean
arrMaxBranchLenSD_crop[i]=getResult("Maximum Branch Length",nResults-3)/scaleFrac; //-- SD
arrLSP_AV_crop[i]=getResult("Longest Shortest Path",nResults-4)/scaleFrac; //-- Mean
arrLSP_SD_crop[i]=getResult("Longest Shortest Path",nResults-3)/scaleFrac; //-- SD
if (verbose>0){print("Number of Results: "+arrCount_crop[i]);}
	} //-- nResults>0
}

//-- Set correct output for per-file data
if (j==0){
outFile=substring(dirList[i], 0, lastIndexOf(dirList[i], "."))+"_alldata.csv";
}
if (j==1){
outFile=substring(dirList[i], 0, lastIndexOf(dirList[i], "."))+"_cropped.csv";
}
saveAs("Results", checkPath(inPath,outFile));
run("Clear Results");

} //-- close main/medial loop

close("*");
} //-- Close file Loop

if (verbose>0){print("------------------------");print("Done processing "+dirList.length+" files");}

//--------------------------------//------------------------------------
// PART 3: Concat and save out the Summary Table
//--------------------------------//------------------------------------
call("ij.plugin.filter.Analyzer.setPrecision", 1);
Array.show("Results (indexes)",dirList,arrCount,arrNumBranchAV,arrNumBranchSD,arrBranchLenAV,arrBranchLenSD,arrMaxBranchLenAV,arrMaxBranchLenSD,arrLSP_AV,arrLSP_SD,arrCount_crop,arrNumBranchAV_crop,arrNumBranchSD_crop,arrBranchLenAV_crop,arrBranchLenSD_crop,arrMaxBranchLenAV_crop,arrMaxBranchLenSD_crop,arrLSP_AV_crop,arrLSP_SD_crop);
saveAs("Results", checkPath(inPath,"Summary.csv"));

restoreSettings(); //-- restore original settings

if (isOpen("Results")) {
	selectWindow("Results");
    run("Close");
}

setBatchMode(false);
print("------------------------");
print("Finished");
print("------------------------");



//-----------------------------------------------------------------------------------------------
// OTHER FUNCTIONS
//-----------------------------------------------------------------------------------------------
function checkPath(path,file){
//-- Check if an output file already exists.
//-- NO: return full path
//-- YES: return full path including a timestamp with grainularity to the second

wholePath=path+File.separator+file;
if (File.exists(wholePath)==0){
	//print("File does not exist");
	return wholePath;
}else{
	//-- Create a unique timestamp (to the second)
	getDateAndTime(year, month, week, day, hour, min, sec, msec);
	//-- Initialse the variable with a string
	ts=d2s(year,0);
	//-- Add the rest on, zero padding - NOTE: MONTH IS ZERO INDEXED
	ts=ts+IJ.pad(month+1,2)+IJ.pad(day,2)+IJ.pad(hour,2)+IJ.pad(min,2)+IJ.pad(sec,2);
	//print("File exists");
	return substring(wholePath, 0, lastIndexOf(wholePath, "."))+"_"+ts+substring(wholePath,lastIndexOf(wholePath, "."));
	}
}