[Threshold HDAB regions using ImageJ] Uses ImageJ's Thresholder to create a selection that gets passed to QuPath as an Annotation or Detection #QuPath #ImageJ #Groovy

Threshold_HDAB_regions_using_ImageJ.groovy

import java.awt.Color
import ij.process.ColorProcessor
import ij.process.AutoThresholder
import ij.plugin.filter.ThresholdToSelection
import ij.IJ

// modified from : 
// https://petebankhead.github.io/qupath/scripting/2018/03/08/script-imagej-to-qupath.html

// This scripts allow threholding of a color deconvoltuion channel ('DAB', 'Hematoxylin',...)

// As the thresholding is done in Imagej (which has an image size limitation), 
// the script works on annotations and we recommend you to process 
// regions having a reasonable size or to use a downsampled version of the image.
// For example, if the pixel size of your image is 0.34 micron, use a value of 1.38  micron 
// for the variable 'requestedPixelSizeMicrons'.
//
// For the threshold, you can either use :
//   - an automatic method that will determine it from the histogram  
//   - a fixed value 
//
// The Annotation(s) to analyse without a name will be named 'Annotation-indexNbr'
//
// The region created from the threshold can be detection or annotation. 
// Using the detection you can make use of the 'Fill detections' button 
// and in combination with the 'Overlay Opacity' slider it will ease viusulisation of the results 
// Measure the area and calculate Area Ratio (Stain/Parent)
//
//
// REQUIRES:
// QuPath 0.4.2
// Last update: Olivier Burri, 20230217


// uncomment to get the ij gui
//IJExtension.getImageJInstance()
IJ.run( "Close All", "" )

// get the micron symbol
um = GeneralTools.micrometerSymbol()

////////////////////////////////////////////////////////////////////////////
//
// PARAMETERS
//
// The stain that should be thresholded
stainName = 'DAB' // 'Hematoxylin' & 'DAB'
stainName = 'Hematoxylin' // 'Hematoxylin' & 'DAB'

// Decicide if scritp uses an Automatic Threholding  Methodd
 useAutoTreshold = false

// Define the method to be used by ImageJ
thresholdMethod = AutoThresholder.Method.Otsu // Default , Moments, Huang, Otsu ...
// ELSE 
// define the FIXED threshold value for the analysis
thresholdMin = 0.1
thresholdMax = 1.0 // color deconvolved images generally have values between 0-1

// Add the classification to the available classifications list
createThClass = true

// Define the resolution of the image we want to work with
// A typical 'full-resolution, 40x' pixel size for a whole slide image is around 0.25 microns
// Therefore we will be requesting at a *much* lower resolution here
requestedPixelSizeMicrons = 1.38

// Sigma value for a Gaussian filter, used to reduce noise before thresholding
// the higher the value the smoother the region
sigmaMicrons = 1 

// Decide if the newly created region should be a detection or an annotation
// The Detections result table will contain the Area of the parent annotation 
// and the ratio DAB-Area / Parent-Area
createDetection = true



////////////////////////////////////////////////////////////////////////////
//
// RUNNER
//



// create  a Class for the thresholded region (allow different coloring than default)
Platform.runLater {
    if ( createThClass ){ 
        if ( stainName == 'DAB' ) createPathClass( stainName, 255,128,0 )
        if ( stainName == 'Hematoxylin' ) createPathClass( stainName, 0,128,255 )
    }
}

// Remove existing regions with stainName
objectList = getAllObjects()

// find the one that contains 'stainName'
def removable = objectList.findAll{ it.getName() =~ stainName }

// remove this list of object and update the Hierarchy
removeObjects( removable , false )


def annotations = getAnnotationObjects()
    
// if there is a selected Annotation, process just that one
selectedObject = getSelectedObject()
if ( selectedObject != null ) annotations = [selectedObject]

annotations.eachWithIndex{ annotation, index -> 

    // Name annotations unless it's already been done
    if ( annotation.getName() == null ){
        annotation.setName( "Annotation-"+index )
    }
    
    process( annotation )
}

fireHierarchyUpdate()

println("Jobs : DONE!")



////////////////////////////////////////////////////////////////////////////
//
// HELPER(s)
//
def createPathClass( def pathClassName, def r, def g, def b ){
    
    // get the List of Class
    available_class = getQuPath().getAvailablePathClasses()
    
    newPathClass = getPathClass( pathClassName )
    
    // check if the new class is already there ,and add it if needed
    if ( !( newPathClass in available_class ) ){
        available_class.add( newPathClass )
    }
    // set the appropriate color
    newPathClass.setColor( getColorRGB( r, g ,b ) )
}



// The main function that does everything
def process( annot ){
    
    // By default, lock the annotation
    annot.setLocked( true )
        
    // Access the relevant QuPath data structures
    def server = getCurrentServer()
    def imageData = getCurrentImageData()
    
    def pixelSize = server.getPixelCalibration().getAveragedPixelSizeMicrons()
    
    // For color deconvolution, we need an 8-bit brightfield RGB image, and also stains to be set
    // Check for these now, and return if we don't have what we need
    def stains = imageData.getColorDeconvolutionStains()
    
    if ( !server.isRGB() || !imageData.isBrightfield() || stains == null ) {
        println 'An 8-bit RGB brightfield image is required!'
        return
    }
    
    // Get the index of the stain we want, based on the specified name
    def selectedStain = stains.getStains( false ).find{ it.getName() == stainName }
    
    def stainIndex = stains.getStainNumber( selectedStain ) - 1
    
    if ( stainIndex < 0 ) {
        println "Could not find stain with name $stainName!"
        return
    }
          
    // Convert requestedPixelSizeMicrons into a sensible downsample value
    int downsample = requestedPixelSizeMicrons / pixelSize
    
    // Create a region request, either for the full image or the selected region
    def region = RegionRequest.createInstance( server.getPath(), downsample, annot.getROI() ) 
    
    // Get the RGB image
    def image = IJTools.convertToImagePlus( server, region )
    def imp = image.getImage()
    //imp.show()
    
    // Get ColorProcessor
    def ip =  imp.getProcessor() as ColorProcessor
    
    def deconvolvedIps = IJTools.colorDeconvolve( ip, stains )

    // Extract the stain ImageProcessor
    def ipStain = deconvolvedIps[ stainIndex ]
   
     // Get the region
    def annot_roi = IJTools.convertToIJRoi( annot.getROI(), image )
    
    // Convert blur sigma to pixels & apply if > 0
    double sigmaPixels = sigmaMicrons / requestedPixelSizeMicrons
    
    if ( sigmaPixels > 0 ) ipStain.blurGaussian( sigmaPixels )
        
    ipStain.setRoi( annot_roi )

    // here we clear outside
    // need to set the color for automatic method on non-rectangle roi
    def color = new Color( 0.0, 0.0, 0.0, 1 ) 
    
    ipStain.setColor( color )
    ipStain.fillOutside( annot_roi )
    
    
    // Set the threshold
    if ( useAutoTreshold ){
        ipStain.setAutoThreshold( thresholdMethod, true )
    } else {
        ipStain.setThreshold( thresholdMin, thresholdMax, 0 )
    }
        
    // Create a selection
    def roiIJ = new ThresholdToSelection().convert( ipStain )
    
    if (roiIJ != null){
        
        ipStain.setRoi(roiIJ)
       
        // Convert ImageJ ROI to a QuPath ROI
        // Here, the pathImage comes in handy because it has the calibration info we want
        def qpROI = IJTools.convertToROI( roiIJ, image )
        
        def thrRegion
        
        // Create a QuPath detection 
        if ( createDetection ) {
            thrRegion = PathObjects.createDetectionObject( qpROI )
        } else {
            thrRegion = PathObjects.createAnnotationObject( qpROI )
        }
        
        thrArea = thrRegion.getROI().getArea() * pixelSize * pixelSize
        annotArea = annot.getROI().getArea() * pixelSize * pixelSize
        areaRatio = thrArea / annotArea
        
        thrRegion.measurements["Area-$stainName $um^2"] = thrArea
        thrRegion.measurements["Area-Parent $um^2"] = annotArea
        thrRegion.measurements["Area Ratio ($stainName / Parent)"] = areaRatio
        
        annot.addChildObject( thrRegion )
    
       // set the name of threshold_region
       def regionName =  stainName
       // after the parent name
       def parentName =  annot.getName().toString()
       if ( parentName != null ) regionName =  parentName+"-"+stainName           
       thrRegion.setName( regionName )
       
       // Optionnal, also set as a Class
       if ( createThClass ) thrRegion.setPathClass( getPathClass( stainName ) )      
    }
}