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


Josef Brandt's avatar
Josef Brandt committed
97
def getParticleStatsWithPixelScale(contour, dataset, fullimage=None, zimg=None, ftir: bool = False):
98 99 100 101 102
    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

Josef Brandt's avatar
Josef Brandt committed
125 126 127 128 129 130 131 132 133 134 135
    if ftir:
        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
Josef Brandt's avatar
Josef Brandt committed
136

137
    return newStats
138

Josef Brandt's avatar
Josef Brandt committed
139

JosefBrandt's avatar
JosefBrandt committed
140 141 142 143 144 145
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
146 147


148 149 150 151 152 153 154 155
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
156

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

Josef Brandt's avatar
Josef Brandt committed
163

164
def getParticleHeight(contour, fullZimg, dataset):
Josef Brandt's avatar
Josef Brandt committed
165
    zimg, _extrema = getParticleImageFromFullimage(contour, fullZimg)
166 167 168 169 170 171 172 173 174 175 176
    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
177

178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193
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
194

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

Josef Brandt's avatar
Josef Brandt committed
199

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

    img = fullimage[ymin:ymax, xmin:xmax]
205
    img = img.copy()
206 207 208 209 210 211 212 213 214 215 216
    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
217 218 219
    return img, (xmin, xmax, ymin, ymax)


220
def contoursToImg(contours, padding=0):
221 222 223 224 225 226
    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))
227 228
    if rangex == 0 or rangey == 0:
        raise InvalidParticleError
229 230 231 232 233 234 235 236 237 238 239

    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
240

241
def imgToCnt(img, xmin, ymin, padding=0):
242
    def getContour(img, contourMode):
243
        if cv2.__version__ > '3.5':
244
            contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
245
        else:
246
            temp, contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
247
        
248 249 250 251 252 253 254
        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)
            
255 256
        return contour
    
257 258 259 260 261 262 263 264 265
    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)
266
    contour = getContour(img, contourMode=cv2.CHAIN_APPROX_NONE)
267
    
268 269 270
    for i in range(len(contour)):
        contour [i][0][0] += xmin-padding
        contour [i][0][1] += ymin-padding
271
        
272
    return contour
273

Josef Brandt's avatar
Josef Brandt committed
274

275 276 277 278 279 280 281
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()
282 283 284 285 286
        ymin, ymax = cnt[:, 1].min(), cnt[:, 1].max()
        
    xmin, xmax = int(round(xmin)), int(round(xmax))
    ymin, ymax = int(round(ymin)), int(round(ymax))
    
287 288
    return xmin, xmax, ymin, ymax

Josef Brandt's avatar
Josef Brandt committed
289

290 291 292 293 294 295 296 297 298
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
299

Josef Brandt's avatar
Josef Brandt committed
300

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

Josef Brandt's avatar
Josef Brandt committed
304

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