particleCharacterization.py 11.1 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, scenePyramid=None, zimg=None, ftir: bool = False):
98 99 100
    if zimg is None:
        zimg = loadZValImageFromDataset(dataset)
        
101
    cnt = deepcopy(contour)
102
    pixelscale = dataset.getPixelScale()
JosefBrandt's avatar
JosefBrandt committed
103
    
104
    newStats = ParticleStats()
105 106 107 108
    newStats.longSize, newStats.shortSize, newStats.area = getContourStats(cnt)
    newStats.longSize *= pixelscale
    newStats.shortSize *= pixelscale
    newStats.area *= (pixelscale**2)
109 110 111 112
    
    if 0 in [newStats.longSize, newStats.shortSize, newStats.area]:
        raise InvalidParticleError
    
113
    newStats.height = getParticleHeight(cnt, zimg, dataset)    
114
    newStats.shape = getParticleShape(cnt, newStats.height)
JosefBrandt's avatar
JosefBrandt committed
115 116 117 118
    if newStats.shape == 'fibre':
        newStats.longSize, newStats.shortSize = getFibreDimension(cnt)
        newStats.longSize *= pixelscale
        newStats.shortSize *= pixelscale
119

Josef Brandt's avatar
Josef Brandt committed
120 121 122 123 124 125
    partImg = None
    if scenePyramid is None and fullimage is not None:
        partImg, extrema = getParticleImageFromFullimage(cnt, fullimage)
    elif scenePyramid is not None and fullimage is None:
        partImg, extrema = getParticleImageFromScenePyramid(cnt, scenePyramid)
    assert partImg is not None, "error in getting particle image"
126
    newStats.color = getParticleColor(partImg)
Josef Brandt's avatar
Josef Brandt committed
127

Josef Brandt's avatar
Josef Brandt committed
128 129 130 131 132 133 134 135 136 137 138
    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
139

140
    return newStats
141

Josef Brandt's avatar
Josef Brandt committed
142

JosefBrandt's avatar
JosefBrandt committed
143 144 145 146 147 148
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
149 150


151 152 153 154 155 156 157 158
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
159

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

Josef Brandt's avatar
Josef Brandt committed
166

167
def getParticleHeight(contour, fullZimg, dataset):
Josef Brandt's avatar
Josef Brandt committed
168
    zimg, _extrema = getParticleImageFromFullimage(contour, fullZimg)
169 170 171 172 173
    if zimg.shape[0] == 0 or zimg.shape[1] == 0:
        raise InvalidParticleError
        
    zimg = cv2.medianBlur(zimg, 5)
    avg_ZValue = np.mean(zimg[zimg > 0])
Josef Brandt's avatar
Josef Brandt committed
174
    if np.isnan(avg_ZValue):  # i.e., only zeros in zimg
175 176 177 178 179
        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
180

181
def getContourStats(cnt):
Josef Brandt's avatar
Josef Brandt committed
182 183 184 185 186 187 188 189 190 191 192 193 194 195
    # 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
196

Josef Brandt's avatar
Josef Brandt committed
197

198 199 200 201
def mergeContours(contours):
    img, xmin, ymin, padding = contoursToImg(contours)
    return imgToCnt(img, xmin, ymin, padding)

Josef Brandt's avatar
Josef Brandt committed
202

203 204 205 206 207
def getParticleImageFromFullimage(contour, fullimage):
    contourCopy = deepcopy(contour)
    xmin, xmax, ymin, ymax = getContourExtrema(contourCopy)

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


Josef Brandt's avatar
Josef Brandt committed
223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242
def getParticleImageFromScenePyramid(contour, scenePyramid):
    contourCopy = deepcopy(contour)
    xmin, xmax, ymin, ymax = getContourExtrema(contourCopy)

    img = scenePyramid.getImagePart(ymin, ymax, xmin, xmax)
    img = img.copy()
    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)
    return img, (xmin, xmax, ymin, ymax)


243
def contoursToImg(contours, padding=0):
244 245
    contourCopy = deepcopy(contours)
    xmin, xmax, ymin, ymax = getContourExtrema(contourCopy)
Josef Brandt's avatar
Josef Brandt committed
246

247 248
    rangex = int(np.round((xmax-xmin)+2*padding))
    rangey = int(np.round((ymax-ymin)+2*padding))
249 250
    if rangex == 0 or rangey == 0:
        raise InvalidParticleError
251 252 253 254 255 256 257 258 259 260 261

    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
262

263
def imgToCnt(img, xmin, ymin, padding=0):
Josef Brandt's avatar
Josef Brandt committed
264

265
    def getContour(img, contourMode):
Josef Brandt's avatar
Josef Brandt committed
266 267 268 269 270 271 272 273 274 275

        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}, '
                  f'taking contour at index {maxIndex}')
            return contours[maxIndex]

276
        if cv2.__version__ > '3.5':
277
            contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
278
        else:
279
            temp, contours, hierarchy = cv2.findContours(img, cv2.RETR_EXTERNAL, contourMode)
280
        
Josef Brandt's avatar
Josef Brandt committed
281
        if len(contours) == 0:  # i.e., no contour found
282
            raise InvalidParticleError
Josef Brandt's avatar
Josef Brandt committed
283
        elif len(contours) == 1:  # i.e., exactly one contour found
284
            contour = contours[0]
Josef Brandt's avatar
Josef Brandt committed
285
        else:  # i.e., multiple contours found
286 287
            contour = getLargestContour(contours)
            
288
        return contour
289 290

    img = closeHolesOfSubImage(img)
291
    contour = getContour(img, contourMode=cv2.CHAIN_APPROX_NONE)
Josef Brandt's avatar
Josef Brandt committed
292

293
    for i in range(len(contour)):
Josef Brandt's avatar
Josef Brandt committed
294 295 296
        contour[i][0][0] += xmin-padding
        contour[i][0][1] += ymin-padding

297
    return contour
298

Josef Brandt's avatar
Josef Brandt committed
299

300 301 302 303 304 305 306
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()
307 308 309 310 311
        ymin, ymax = cnt[:, 1].min(), cnt[:, 1].max()
        
    xmin, xmax = int(round(xmin)), int(round(xmax))
    ymin, ymax = int(round(ymin)), int(round(ymax))
    
312 313
    return xmin, xmax, ymin, ymax

Josef Brandt's avatar
Josef Brandt committed
314

315 316 317 318 319 320 321 322 323
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
324

Josef Brandt's avatar
Josef Brandt committed
325

326
def loadFullimageFromDataset(dataset):
327
    return cv2imread_fix(dataset.getImageName())
328

Josef Brandt's avatar
Josef Brandt committed
329

330 331
def loadZValImageFromDataset(dataset):
    return cv2imread_fix(dataset.getZvalImageName(), cv2.IMREAD_GRAYSCALE)