// Optimized hiPSC Co Culture in vitro segmentation - Bardy Lab - SAHMRI

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

// To run this script, click 'Run' below and you will be prompted to load the first image in your experiments image sequence. 
// 		1) Double-click on JUST THE FIRST IMAGE in your experiment
// 		2) Check off the "sort numerically" and 
// 		2) "load virtual stack" options, 
// 		3) click OK, 
// and the script will automatically do the rest. (perform some processing steps and generate a list of ROIs)

// This script will take some time to process. The most time consuming part at this point is recording the ROIs, which is limited by data transfer speed. 
// I wrote the script so that you can run it on a low-RAM machine if need be. // But if your files are being processed on a local machine (not through a server) then it should read much faster. 

// The highlighted ROIs should be most of the cells in the field of view.
// Somewhat unbiased to whether the cells are responding or not. 

// Several items will be saved and output into each experiments folder
// 		1) Cell Timeseries: output as a .csv file in \results folder
// 				Each collumn = cell/roi, each row = timepoint
// 				May be too large to open in Excel. 
// 				if you need spatial locations, stdev,etc  its easy to add
// 		2) ROIs: this is saved as a .zip file in the \results folder
// 		3) AvgIP - an averaged image for segmenting saved in \zproj

// I attempted to comment the script as throroughly as possible so you can figure out whats going on at each step of the way - tweaking it as needed for different application. Feel free to mess around with it

// In publications, the papers to cite using for this scripts are at the bottom of this page. Those papers are for FIJI and MorphoLibJ.

// Always happy to help explain how the processing works as needed, as well as help fix any issues you run into. 

// ==================================================================
// Written be Wilson R Adams | Vanderbilt Biophotonics Center | Dec 2019
// ==================================================================

// For Bardy lab, load the first raw image of your timeseries stack with this command:
// ================================================
roiManager("reset");
run("Image Sequence...");
// ================================================

// For Wilson, Load images with this command:
// ================================================
// run("TIFF Virtual Stack...");
// =================================================


orig = getTitle();
filename = File.nameWithoutExtension();	
dir = getDirectory("image");
filelist = getFileList(dir);
Array.show(filelist);

// Make all of the save directories for results
zpd = dir+"\\zproj";
File.makeDirectory(zpd);
crp = dir+"\\crop"; 
File.makeDirectory(crp);
resdir = dir+"\\results"; 
File.makeDirectory(resdir);

imgc = getTitle();
run("Enhance Contrast", "saturated=0.5");

// Average Z-projections to improve SNR for segmentation
run("Z Project...", "start=5 stop=60 projection=[Average Intensity]");
// Crop out heading
makeRectangle(0, 9, 2048, 2040);
run("Crop");
avgip = getTitle();	// ** SAVE AVGIP **

// get stack dimensions
Stack.getDimensions(width, height, channels, slices, frames);

// Flatten out image field with background subtraction and CLAHE
run("Subtract Background...", "rolling=50");
run("Enhance Contrast", "saturated=1");
run("Green Fire Blue ");
run("Enhance Local Contrast (CLAHE)", "blocksize=100 histogram=256 maximum=20.00 mask=*None* fast_(less_accurate)");

// Enhance cell bodies and suppress dendtritic projection features with adisk-shaped morphological opening filter
selectWindow(avgip)
run("Morphological Filters", "operation=Opening element=Disk radius=8");
run("Gaussian Blur...", "sigma=4");
rename("filt");
filt = getTitle();

// Segment out optical fiber (if you need to... otherwise comment out). 
//waitForUser("Use the Rectangle ROI to cut out the fiber from the FOV")

// Locate cell bodies as local maxima and record
run("Find Maxima...", "prominence=10 output=[Point Selection]");
roiManager("add");

// Make Marker image for watershed with local maxima points
newImage("markers", "8-bit black", width, height, 1);
markers = getTitle();
roiManager("Select", 0);
roiManager("Fill");
run("Make Binary");

// Threshold filtered iamge to get mask for watershed mask
selectWindow(filt)
run("Duplicate...", " ");
setAutoThreshold("Huang dark");
run("Threshold...");
run("Convert to Mask");
rename("mask");
mask = getTitle();  // ** SAVE MASK **

// Get gradient image of mask for watershed input
run("Morphological Filters", "operation=Gradient element=Disk radius=2");
rename("grad");
grad = getTitle();

// Run watershed Segmentation
run("Marker-controlled Watershed", "input=grad marker=markers mask=mask binary calculate use");
wtsh = getTitle(); 
run("Fire");
roiManager("reset");

run("Duplicate...", " ");
setThreshold(0.1000, 1000000000000000000000000000000.0000);
run("Convert to Mask");
wsmsk = getTitle();

// Circle cell data
run("Analyze Particles...", "size=250-Infinity pixel circularity=0.10-1.00 add");
close(wsmsk);

// RUN CODE ABOVE THIS LINE TO GET ROIS WITHOUT MEASURING AND SAVING
// 		[	highlight code and press Ctrl+Shift+R	]
// =================================================
// =================================================
// =================================================
// =================================================

// Apply ROIs to Original image
selectWindow(orig)
roiManager("Show None");
roiManager("Show All");

// Measure the timeseries from original raw data stack
roiManager("multi-measure one");	// **SAVE RESULTS **

// Save Results and ROIs
resname = "results_"+orig;
print(resname);
restitle = filename;
saveAs("Results", resdir+'\\'+restitle+"_results.csv");

roiManager("Save", resdir+'\\'+restitle+"_roi.zip");
close("Results");
roiManager("reset");

// Save AvgIP image
selectWindow(avgip)
rename("avgip_"+filename);
roiManager("Show None");
roiManager("Show All");
saveAs("tiff", zpd+"\\"+getTitle());
avgip = getTitle();


// CLOSE stuff
close(markers);
close(grad);
close(wtsh);
close(mask);
close(filt)
//close(orig);
//close(avgip);


// PUBLICATIONS: 
// FIJI Publication:
// 	Schindelin, J.; Arganda-Carreras, I. & Frise, E. et al. (2012), "Fiji: an open-source platform for biological-image analysis", Nature methods 9(7): 676-682, PMID 22743772, doi:10.1038/nmeth.2019 

// MorphoLibJ Publication:
// Legland, D.; Arganda-Carreras, I. & Andrey, P. (2016), "MorphoLibJ: integrated library and plugins for mathematical morphology with ImageJ", Bioinformatics (Oxford Univ Press) 32(22): 3532-3534, PMID 27412086, doi:10.1093/bioinformatics/btw413


// Interactive data viewing


// Make a pretty video
// AvgIP --> Create Background --> Subtract from denoised stack
//