particleCharacterization.py 10.3 KB
Newer Older
JosefBrandt's avatar
JosefBrandt committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
GEPARD - Gepard-Enabled PARticle Detection
Copyright (C) 2018  Lars Bittrich and Josef Brandt, Leibniz-Institut für 
Polymerforschung Dresden e. V. <bittrich-lars@ipfdd.de>    

This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License
along with this program, see COPYING.  
If not, see <https://www.gnu.org/licenses/>.
"""

import numpy as np
import cv2
25
from copy import deepcopy
JosefBrandt's avatar
JosefBrandt committed
26

27 28
from .particleClassification.colorClassification import ColorClassifier
from .particleClassification.shapeClassification import ShapeClassifier
Josef Brandt's avatar
Josef Brandt committed
29
from .ftirAperture import findRotatedMaximalRectangle
30 31
from ..segmentation import closeHolesOfSubImage
from ..errors import InvalidParticleError
32
from ..helperfunctions import cv2imread_fix
33

Josef Brandt's avatar
Josef Brandt committed
34 35 36 37 38 39 40 41 42 43 44 45 46 47

class FTIRAperture(object):
    """
    Configuration for an FTIR aperture. CenterCoords, width and height in Pixels
    """
    centerX: float = None
    centerY: float = None
    centerZ: float = None
    width: float = None
    height: float = None
    angle: float = None
    rectCoords: list = None


48
class ParticleStats(object):
Josef Brandt's avatar
Josef Brandt committed
49 50 51 52 53 54 55 56
    longSize: float = None
    shortSize: float = None
    height: float = None
    area: float = None
    shape: str = None
    color: str = None
    aperture: FTIRAperture = None

57 58 59 60 61 62 63 64

def particleIsValid(particle):
    if particle.longSize == 0 or particle.shortSize == 0:
        return False
    
    if cv2.contourArea(particle.contour) == 0:
        return False
    return True
Josef Brandt's avatar
Josef Brandt committed
65 66


67 68 69 70 71 72 73 74 75 76 77 78 79
def updateStatsOfParticlesIfNotManuallyEdited(particleContainer):
    dataset = particleContainer.datasetParent
    fullimage = loadFullimageFromDataset(dataset)
    zimg = loadZValImageFromDataset(dataset)
    
    for particle in particleContainer.particles:
        if not hasattr(particle, "wasManuallyEdited"):
            particle.wasManuallyEdited = False
        
        if not particle.wasManuallyEdited:
            newStats = getParticleStatsWithPixelScale(particle.contour, dataset, fullimage, zimg)
            particle.__dict__.update(newStats.__dict__)

Josef Brandt's avatar
Josef Brandt committed
80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96

def getFTIRAperture(partImg: np.ndarray) -> FTIRAperture:
    rectPoints, angle, width, heigth = findRotatedMaximalRectangle(partImg, nbre_angle=4, limit_image_size=300)
    xvalues: list = [i[0] for i in rectPoints]
    yvalues: list = [i[1] for i in rectPoints]

    aperture = FTIRAperture()
    aperture.height = heigth
    aperture.width = width
    aperture.angle = angle
    aperture.centerX = int(round(min(xvalues) + (max(xvalues) - min(xvalues)) / 2))
    aperture.centerY = int(round(min(yvalues) + (max(yvalues) - min(yvalues)) / 2))
    aperture.rectCoords = rectPoints

    return aperture


97 98 99 100 101 102
def getParticleStatsWithPixelScale(contour, dataset, fullimage=None, zimg=None):
    if fullimage is None:
        fullimage = loadFullimageFromDataset(dataset)
    if zimg is None:
        zimg = loadZValImageFromDataset(dataset)
        
103
    cnt = deepcopy(contour)
104
    pixelscale = dataset.getPixelScale()
JosefBrandt's avatar
JosefBrandt committed
105
    
106
    newStats = ParticleStats()
107 108 109 110
    newStats.longSize, newStats.shortSize, newStats.area = getContourStats(cnt)
    newStats.longSize *= pixelscale
    newStats.shortSize *= pixelscale
    newStats.area *= (pixelscale**2)
111 112 113 114
    
    if 0 in [newStats.longSize, newStats.shortSize, newStats.area]:
        raise InvalidParticleError
    
115
    newStats.height = getParticleHeight(cnt, zimg, dataset)    
116
    newStats.shape = getParticleShape(cnt, newStats.height)
JosefBrandt's avatar
JosefBrandt committed
117 118 119 120
    if newStats.shape == 'fibre':
        newStats.longSize, newStats.shortSize = getFibreDimension(cnt)
        newStats.longSize *= pixelscale
        newStats.shortSize *= pixelscale
121

Josef Brandt's avatar
Josef Brandt committed
122
    partImg, extrema = getParticleImageFromFullimage(cnt, fullimage)
123
    newStats.color = getParticleColor(partImg)
Josef Brandt's avatar
Josef Brandt committed
124 125 126 127 128 129 130 131 132 133 134 135

    padding: int = int(2)
    imgforAperture = np.zeros((partImg.shape[0]+2*padding, partImg.shape[1]+2*padding))
    imgforAperture[padding:partImg.shape[0]+padding, padding:partImg.shape[1]+padding] = np.mean(partImg, axis=2)
    imgforAperture[imgforAperture > 0] = 1  # convert to binary image
    imgforAperture = np.uint8(1 - imgforAperture)  # background has to be 1, particle has to be 0
    aperture: FTIRAperture = getFTIRAperture(imgforAperture)
    aperture.centerX = aperture.centerX + extrema[0] - padding
    aperture.centerY = aperture.centerY + extrema[2] - padding
    aperture.centerZ = newStats.height
    newStats.aperture = aperture

136
    return newStats
137

Josef Brandt's avatar
Josef Brandt committed
138

JosefBrandt's avatar
JosefBrandt committed
139 140 141 142 143 144
def getFibreDimension(contour):
    longSize = cv2.arcLength(contour, True)/2
    img = contoursToImg([contour])[0]
    dist = cv2.distanceTransform(img, cv2.DIST_L2, 3)
    maxThickness = np.max(dist)*2
    return longSize, maxThickness
Josef Brandt's avatar
Josef Brandt committed
145 146


147 148 149 150 151 152 153 154
def getParticleColor(imgRGB, colorClassifier=None):
    img = cv2.cvtColor(imgRGB, cv2.COLOR_RGB2HSV_FULL)
    meanHSV = cv2.mean(img)
    if colorClassifier is None:
        colorClassifier = ColorClassifier()
    color = colorClassifier.classifyColor(meanHSV)
    return color

Josef Brandt's avatar
Josef Brandt committed
155

156 157 158
def getParticleShape(contour, particleHeight, shapeClassifier=None):
    if shapeClassifier is None:
        shapeClassifier = ShapeClassifier()
JosefBrandt's avatar
JosefBrandt committed
159
    shape = shapeClassifier.classifyShape(contour, particleHeight)
160 161
    return shape

Josef Brandt's avatar
Josef Brandt committed
162

163
def getParticleHeight(contour, fullZimg, dataset):
Josef Brandt's avatar
Josef Brandt committed
164
    zimg, _extrema = getParticleImageFromFullimage(contour, fullZimg)
165 166 167 168 169 170 171 172 173 174 175
    if zimg.shape[0] == 0 or zimg.shape[1] == 0:
        raise InvalidParticleError
        
    zimg = cv2.medianBlur(zimg, 5)
    avg_ZValue = np.mean(zimg[zimg > 0])
    if np.isnan(avg_ZValue):  #i.e., only zeros in zimg
        avg_ZValue = 0
    z0, z1 = dataset.zpositions.min(), dataset.zpositions.max()
    height = avg_ZValue/255.*(z1-z0) + z0
    return height

Josef Brandt's avatar
Josef Brandt committed
176

177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192
def getContourStats(cnt):
    ##characterize particle
        if cnt.shape[0] >= 5:       ##at least 5 points required for ellipse fitting...
            ellipse = cv2.fitEllipse(cnt)
            short, long = ellipse[1]
        else:
            rect = cv2.minAreaRect(cnt)
            long, short = rect[1]
            
        if short>long:
            long, short = short, long
        
        area = cv2.contourArea(cnt)
    
        return long, short, area

Josef Brandt's avatar
Josef Brandt committed
193

194 195 196 197
def mergeContours(contours):
    img, xmin, ymin, padding = contoursToImg(contours)
    return imgToCnt(img, xmin, ymin, padding)

Josef Brandt's avatar
Josef Brandt committed
198

199 200 201 202 203
def getParticleImageFromFullimage(contour, fullimage):
    contourCopy = deepcopy(contour)
    xmin, xmax, ymin, ymax = getContourExtrema(contourCopy)

    img = fullimage[ymin:ymax, xmin:xmax]
204
    img = img.copy()
205 206 207 208 209 210 211 212 213 214 215
    mask = np.zeros(img.shape[:2])
    
    for i in range(len(contourCopy)):
        contourCopy[i][0][0] -= xmin
        contourCopy[i][0][1] -= ymin
        
    cv2.drawContours(mask, [contourCopy], -1, (255, 255, 255), -1)
    cv2.drawContours(mask, [contourCopy], -1, (255, 255, 255), 1)
    
    img[mask == 0] = 0
    img = np.array(img, dtype = np.uint8)
Josef Brandt's avatar
Josef Brandt committed
216 217 218
    return img, (xmin, xmax, ymin, ymax)


219
def contoursToImg(contours, padding=0):
220 221 222 223 224 225
    contourCopy = deepcopy(contours)
    xmin, xmax, ymin, ymax = getContourExtrema(contourCopy)
    
    padding = padding   #pixel in each direction
    rangex = int(np.round((xmax-xmin)+2*padding))
    rangey = int(np.round((ymax-ymin)+2*padding))
226 227
    if rangex == 0 or rangey == 0:
        raise InvalidParticleError
228 229 230 231 232 233 234 235 236 237 238

    img = np.zeros((rangey, rangex))
    for curCnt in contourCopy:
        for i in range(len(curCnt)):
            curCnt[i][0][0] -= xmin-padding
            curCnt[i][0][1] -= ymin-padding
        cv2.drawContours(img, [curCnt], -1, 255, -1)
        cv2.drawContours(img, [curCnt], -1, 255, 1)
    img = np.uint8(cv2.morphologyEx(img, cv2.MORPH_CLOSE, np.ones((3, 3))))
    return img, xmin, ymin, padding

Josef Brandt's avatar
Josef Brandt committed
239

240
def imgToCnt(img, xmin, ymin, padding=0):
241
    def getContour(img, contourMode):
242
        if cv2.__version__ > '3.5':
243
            contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
244
        else:
245
            temp, contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
246
        
247 248 249 250 251 252 253
        if len(contours) == 0:   #i.e., no contour found
            raise InvalidParticleError
        elif len(contours) == 1:   #i.e., exactly one contour found
            contour = contours[0]
        else:       #i.e., multiple contours found
            contour = getLargestContour(contours)
            
254 255
        return contour
    
256 257 258 259 260 261 262 263 264
    def getLargestContour(contours):
        areas = []
        for contour in contours:
            areas.append(cv2.contourArea(contour))
        maxIndex = areas.index(max(areas))
        print(f'{len(contours)} contours found, getting the largest one. Areas are: {areas}, taking contour at index {maxIndex}')
        return contours[maxIndex]

    img = closeHolesOfSubImage(img)
265
    contour = getContour(img, contourMode=cv2.CHAIN_APPROX_NONE)
266
    
267 268 269
    for i in range(len(contour)):
        contour [i][0][0] += xmin-padding
        contour [i][0][1] += ymin-padding
270
        
271
    return contour
272

Josef Brandt's avatar
Josef Brandt committed
273

274 275 276 277 278 279 280
def getContourExtrema(contours):
    try:
        cnt = np.vstack(tuple(contours))
        xmin, xmax = cnt[:,0,:][:, 0].min(), cnt[:,0,:][:, 0].max()
        ymin, ymax = cnt[:,0,:][:, 1].min(), cnt[:,0,:][:, 1].max()        
    except IndexError: #i.e., not a list of contours was passed, but an individual contour. Hence, the above indexing does not work
        xmin, xmax = cnt[:, 0].min(), cnt[:, 0].max()
281 282 283 284 285
        ymin, ymax = cnt[:, 1].min(), cnt[:, 1].max()
        
    xmin, xmax = int(round(xmin)), int(round(xmax))
    ymin, ymax = int(round(ymin)), int(round(ymax))
    
286 287
    return xmin, xmax, ymin, ymax

Josef Brandt's avatar
Josef Brandt committed
288

289 290 291 292 293 294 295 296 297
def getParticleCenterPoint(contour):
    img, xmin, ymin, padding = contoursToImg(contour)
    dist = cv2.distanceTransform(img, cv2.DIST_L2, 3)
    ind = np.argmax(dist)
    y = ind//dist.shape[1]-1
    x = ind%dist.shape[1]-1
    x += xmin
    y += ymin
    return x, y
298

Josef Brandt's avatar
Josef Brandt committed
299

300
def loadFullimageFromDataset(dataset):
301
    return cv2imread_fix(dataset.getImageName())
302

Josef Brandt's avatar
Josef Brandt committed
303

304 305
def loadZValImageFromDataset(dataset):
    return cv2imread_fix(dataset.getZvalImageName(), cv2.IMREAD_GRAYSCALE)