// MANTIS_AutoExtract_multimodeprocess.ijm

// MANTIS Auto-extract Macro with Multimode Data Processing 
// and Extraction

// Automatically imports, names, and saves condensed image data from MANTIS into the parent directory for subsequent processing. 

// Processing pipeline from multimode_seg_extract.ijm is applied

// Process is repeated for all folders in the directory


// Start by navigating to the parent directory of your experiments are stored when prompted and let the macro run. 

// Written by Wilson R Adams, 23 July 2020
 // ======================================================

// Get file directory from user input
dir = getDirectory("Navigate to the folder to be processed");
print("STARTING file extraction: "+dir);

filelist = getFileList(dir);
Array.show(filelist);

// Setup Workspace
roiManager("reset");

// Loop through each folder in the parent directory

//for (i = 0; i < 1; i++) {			// for debugging one folder at a time
for (i = 0; i < filelist.length; i++) {
    // only open folders, not files in the parent direcory
    if (endsWith(filelist[i], "/")) {  
    	
		// List image files in sub folder
		sublist = getFileList(filelist[i]);				 
		firsttif = dir+filelist[i]+"ChanA_0001_0001_0001_0001.tif";
		
		// If folder has images in it, 
		if (File.exists(firsttif)) {
			// Import image using Bioformats to preserve metadata. 
			print("Loading images...");
			
		// Uses virtual stacks for speed and RAM, especially on server
			run("Bio-Formats Importer", "open=[&firsttif] color_mode=Default rois_import=[ROI manager] split_channels view=Hyperstack stack_order=XYCZT use_virtual_stack");

			// Rename all of the files once imported (assumes MANTIS)
			print("Renaming Images Appropriately");
			runMacro("MANTIS_imgrename_ext.ijm");
			
			// Save in parent directory and close all of the images
			wins = getList("image.titles");
			print("Saving Images....");
			
			for (j = 0; j < wins.length; j++) {
				selectWindow(wins[j]);
				save(dir+wins[j]+".tif");
			}
			// Run the Auto Segmentation and Data Extraction 
			print("Starting Segmentation and Data Extraction...");
			multimode_seg_extract(wins[1], wins[0], dir);
			close("*");
			// Print progress once images are extracted
			print(filelist[i]+" COMPLETE");
		} // end of 'file exist' if stmt

    } // end of tile folder if statment
   
} // end of image folder for loop
print("FINISHED file extraction: "+dir);


function multimode_seg_extract(seg, func, dir) { 
// This function segments the stack 'seg' and applies the ROIs to 
// Functional image stack 'func' and 'seg' while saving results to 
// directory 'dir'

// Segmentation done with modified AstroSeg_Virtual routine

// All data saved to be imported and analyzed in program of choice

// Written by Wilson Adams, Vanderbilt Biophotonics Center, July 2020

srs = seg; 	// Segmentable Image Stack
via = func;	// Functional Image Stack (Viability, Calcium, etc)

// Get ROIs from AvgIP of SRS stack
multimode_seg_extracter(srs);
srsavgip = getTitle();

// Get Centroids of ROIs for location transform
run("Set Measurements...", "centroid redirect=None decimal=6");
roiManager("multi measure one");
saveAs("Results", dir+srs+"_centroid");
roiManager("save", dir+srs+"_roi.zip");

// Get timeseries
run("Set Measurements...", "mean redirect=None decimal=6");

vianame = File.getName(via);
selectWindow(via);
roiManager("show all");
roiManager("multi-measure one");
saveAs("Results", dir+vianame+"_timeseries");
close("Results");

srsname = File.getName(srs);
selectWindow(srsname);
roiManager("show all");
roiManager("multi-measure one");
saveAs("Results", dir+srsname+"_timeseries");
close("Results");

selectWindow(srsavgip);
saveAs("tif", dir+srsname+"_avgip"+".tif");

roiManager("reset");
}

// ============================================

function multimode_seg_extracter(segstack) { 
	// given a segmentable (segstack) and functional (funcstack) 
	// image stacks,this macro will segment cells from image 

	// step 1 - get segmentable image (Avg IP from seg stack)
	selectWindow(segstack);
	run("Z Project...", "projection=[Average Intensity]");
	segAvgIP = getTitle();

	// Step 2 - run segmenting code on it to get ROIs. 
	AstroSeg_virtual(segAvgIP); // Modify function below
}




// ===========================================

// Run on AvgIP to debug
//test = getTitle();
//AstroSeg_virtual(test);

function AstroSeg_virtual(avgip) { 
// Modified from Astrocyte Segmentation

// Open AvgIP file, get image stats
selectWindow(avgip);
Stack.getDimensions(width, height, channels, slices, frames);
getPixelSize(unit, pixelWidth, pixelHeight);

// Sub Background (50px kernel)
run("Duplicate...", " ");
wkimg = getTitle();
// Gaussian Blur (7px kernel)
run("Gaussian Blur...", "sigma=2");
run("Median...", "radius=1");
run("Subtract Background...", "rolling=30");
run("Enhance Contrast", "saturated=5.00");




// Additional Flatfield Correction if needed for cell activation
run("Enhance Local Contrast (CLAHE)", "blocksize=20 histogram=256 maximum=12.00 mask=*None* fast_(less_accurate)");
run("Duplicate...", " ");

// =============== USER MODIFY ===============
// Threshold image to generate Mask (Otsu(astro) - Li (BV2)). Erode a bit. Fill holes.
setAutoThreshold("Huang dark");
//run("Threshold...");
run("Convert to Mask");
rename("mask");
mask = getTitle();

// Make Gradient mask image
run("Morphological Filters", "operation=Gradient element=Disk radius=2");
rename("grad");
grad = getTitle();

// =============== USER MODIFY ===============
// Find Local Maxima (Prominence 20 - Astro | 55 BV2), add to ROI manager
selectWindow(wkimg);
run("Find Maxima...", "prominence=10 output=[Point Selection]");
roiManager("add");

// Make blank marker image. Apply points to img with 'Fill'
// You cant use ROImanager point selections. you have to have an image with points identifies as maximum seed points.
newImage("markers", "8-bit black", width, height, 1);
markers = getTitle();
roiManager("Select", 0);
roiManager("Fill");
run("Make Binary");
run("Properties...", "channels=1 slices=1 frames=1 unit=&unit pixel_width=&pixelWidth pixel_height=&pixelHeight voxel_depth=1.0000");

// Marker based watershed (input = gradient, binary markers, mask)
run("Marker-controlled Watershed", "input=grad marker=markers mask=mask binary calculate use");
wtsh = getTitle(); 
run("Fire");
roiManager("select", 0);
roiManager("delete");

// Threshold watershed map 1:inf for Analyze Particles
run("Duplicate...", " ");
setThreshold(0.1000, 1000000000000000000000000000000.0000);
run("Convert to Mask");
wsmsk = getTitle(); 

// Analyze particles --> min area = 200 into ROI. Save
run("Analyze Particles...", "size=10-Infinity pixel circularity=0.10-1.00 add");

// Check and Update ROIs (if needed, probably needed...). Save to new subdir

selectWindow(avgip); 
roiManager("show all with labels");

close(wkimg);
close(mask);
close(grad);
close(markers);
close(wtsh);
close(wsmsk);
}