[Crop Dataset Around Nuclei] Crops an Imaris Dataset where nuclei were tracked by taking the centroid of each surface at each timepoint #groovy #imaris #cool

Crop Dataset Around Nuclei.groovy

/** 
 * Split Imaris Tracked Nuclei into separate ims files 
 * 
 * 
 * Author: Olivier Burri EPFL SV PTECH BIOP
 * for Noémie Chabot, Vastenhouw Lab
 * 
 * Last update: 01-10-2021
 * 
 * Dependencies: EasyXT https://github.com/BIOP
 * 
 * Inputs
 * Imaris with open dataset that has tracked surfaces called "Nucleus"
 * 
 * Outputs:
 * One Imaris file per track centered around the centroid of each surface for each timepoint
 * 
 */

#@ File outputFolder (label="Output Folder", style="directory")
#@ Double extentXYMicrons (label="Size in XY for cropped volumes (um)", value=15)


// Imports to use the different plugins
import ij.*
import ch.epfl.biop.imaris.*
import ij.plugin.Duplicator
import ij.plugin.*
import Imaris.*;
import Imaris.tType;

// Cleanup for ImageJ
IJ.run("Close All", "")
EasyXT.Utils.resetImarisConnection()

// Get the name of the currently open image
def name = EasyXT.Files.getOpenFileName()

println( "Processing Image "+name)

// Make a copy of the current dataset from which we will crop out the different nuclei
def raw_data =  EasyXT.Dataset.getCurrent().Clone()

// Get a reference to the tracked nuclei surfaces previously created by the user
def surfaces = EasyXT.Surfaces.find("Nucleus")


// Get the track information, namely the track IDs and the associated surface IDs for each track
def edges = surfaces.GetTrackEdges()
def trackIDs = surfaces.GetTrackIds()

// Composite them in order to get a meaningful list with TrackID linked to each surface ID edge
def eid = [edges, trackIDs].transpose().collect{ e, id -> 
	return [id, e]
}

// Find how many tracks we have by keeping only the unique TrackIDs
uniqueIds = (trackIDs as List).toUnique()

// Get Calibration for knowing how to convert between um and px
def cal = EasyXT.Dataset.getCalibration()

// Define the size of the cropped image volumes
w = cal.getRawX(extentXYMicrons) as int
h = cal.getRawY(extentXYMicrons) as int
d = raw_data.GetSizeZ()
nC = raw_data.GetSizeC()

// Remove the surfaces as otherwise they will show up on each image, which is ugly
EasyXT.Scene.removeItem(surfaces)

// Loop all TrackIds and get all of the IDs of the surfaces from the edge map
uniqueIds.each{ tID ->

	// Find all surface IDs linked to this particular track and keep only the unique ones
	def surfaceIds = eid.findAll{ it[0].equals(tID) }.collect{ it[1] }.flatten().toUnique()

	// Prepare the dataset that will contain the cropped nucleus track
	// This is ugly but EasyXT does not support it yet
	def cropped = EasyXT.Utils.getImarisApp().GetFactory().CreateDataSet()
	cropped.Create(tType.eTypeUInt16 , w, h, d, 2, raw_data.GetSizeT())
	
	// Make sure that the channel names and colors are ok
	(0..nC-1).each{ c ->
		cropped.SetChannelColorRGBA(c, raw_data.GetChannelColorRGBA(c))
		cropped.SetChannelName(c, raw_data.GetChannelName(c))
		cropped.SetChannelRange (c, raw_data.GetChannelRangeMin (c), raw_data.GetChannelRangeMax(c))
	}

	// More importantly set the calibration in imaris by setting the exents and the time separation	
	cropped.SetExtendMinX(0)
	cropped.SetExtendMinY(0)
	cropped.SetExtendMinZ(0)
	cropped.SetExtendMaxX(cal.pixelWidth * w as float)
	cropped.SetExtendMaxY (cal.pixelWidth * h as float)
	cropped.SetExtendMaxZ (cal.pixelDepth *d as float)
	cropped.SetTimePointsDelta( raw_data.GetTimePointsDelta() )

	
	// Loop through ids , get center coordinates
	surfaceIds.each{ sid ->

		// Something easy!
		def cm = surfaces.GetCenterOfMass(sid).flatten()

		
		def t = surfaces.GetTimeIndex(sid)

		// get the centroid as pixel coordinates for the extraction
		def x = cm[0] / cal.pixelWidth as int
		def y = cm[1] / cal.pixelWidth as int
		def z = cm[2] / cal.pixelDepth as int
	
		// Create box around data for cropping
		def xpx = x-w/2 as int
		def ypx = y-h/2 as int
		
		// Loop through the channels for cropping
		(0..nC-1).each{ c ->
			// Ge the volume
			def data_t = raw_data.GetDataSubVolumeShorts(xpx, ypx, 0, c, t, w, h, d)

			//set the volume
			cropped.SetDataSubVolumeShorts( data_t, 0, 0, 0,c ,t )
		}

		//Make sure that the timepoint is correct otherwise the metadata is wrong
		cropped.SetTimePoint (t, raw_data.GetTimePoint(t))
		
	} // end of surfaceID loop

	// set the dataset to the Iamris Scene for saving
	EasyXT.Dataset.setCurrent(cropped)

	// Save using the trackID
	def savedFile = new File( outputFolder, ""+tID+name )
	EasyXT.Files.saveImage(savedFile)
}

// Mop up: make current dataset as it was before
EasyXT.Dataset.setCurrent(raw_data)
//EasyXT.Scene.addIdtem(surfaces)


println ("Done")