MaximumCurvature.py 8.3 KB
Newer Older
Pedro TOME's avatar
Pedro TOME committed
1 2 3
#!/usr/bin/env python
# vim: set fileencoding=utf-8 :

4 5 6
import math
import numpy

Pedro TOME's avatar
Pedro TOME committed
7 8 9
import bob.core
import bob.io.base

10
from bob.bio.base.extractor import Extractor
11

Pedro TOME's avatar
Pedro TOME committed
12 13
from .. import utils

14

Pedro TOME's avatar
Pedro TOME committed
15
class MaximumCurvature (Extractor):
16 17 18 19 20 21
  """MiuraMax feature extractor

  Based on N. Miura, A. Nagasaka, and T. Miyatake, Extraction of Finger-Vein
  Pattern Using Maximum Curvature Points in Image Profiles. Proceedings on IAPR
  conference on machine vision applications, 9 (2005), pp. 347--350

22 23 24
  Parameters:

    sigma (int, Optional): Sigma used for determining derivatives
25

26
  """
Pedro TOME's avatar
Pedro TOME committed
27

28

29 30
  def __init__(self, sigma = 5):
    Extractor.__init__(self, sigma = sigma)
Pedro TOME's avatar
Pedro TOME committed
31
    self.sigma = sigma
32 33 34 35 36 37


  def maximum_curvature(self, image, mask):
    """Computes and returns the Maximum Curvature features for the given input
    fingervein image"""

Pedro TOME's avatar
Pedro TOME committed
38 39 40 41 42
    if image.dtype != numpy.uint8:
       image = bob.core.convert(image,numpy.uint8,(0,255),(0,1))
    #No es necesario pasarlo a uint8, en matlab lo dejan en float64. Comprobar si varian los resultados en vera database y ajustar.

    finger_mask = numpy.zeros(mask.shape)
43 44
    finger_mask[mask == True] = 1

Pedro TOME's avatar
Pedro TOME committed
45
    winsize = numpy.ceil(4*self.sigma)
46

Pedro TOME's avatar
Pedro TOME committed
47 48 49
    x = numpy.arange(-winsize, winsize+1)
    y = numpy.arange(-winsize, winsize+1)
    X, Y = numpy.meshgrid(x, y)
50

Pedro TOME's avatar
Pedro TOME committed
51 52 53 54 55 56
    h = (1/(2*math.pi*self.sigma**2))*numpy.exp(-(X**2 + Y**2)/(2*self.sigma**2))
    hx = (-X/(self.sigma**2))*h
    hxx = ((X**2 - self.sigma**2)/(self.sigma**4))*h
    hy = hx.T
    hyy = hxx.T
    hxy = ((X*Y)/(self.sigma**4))*h
57

Pedro TOME's avatar
Pedro TOME committed
58
    # Do the actual filtering
59

60 61 62 63 64
    fx = utils.imfilter(image, hx)
    fxx = utils.imfilter(image, hxx)
    fy = utils.imfilter(image, hy)
    fyy = utils.imfilter(image, hyy)
    fxy = utils.imfilter(image, hxy)
65

Pedro TOME's avatar
Pedro TOME committed
66 67 68 69
    f1  = 0.5*numpy.sqrt(2)*(fx + fy)   # \  #
    f2  = 0.5*numpy.sqrt(2)*(fx - fy)   # /  #
    f11 = 0.5*fxx + fxy + 0.5*fyy       # \\ #
    f22 = 0.5*fxx - fxy + 0.5*fyy       # // #
70

Pedro TOME's avatar
Pedro TOME committed
71
    img_h, img_w = image.shape  #Image height and width
72

Pedro TOME's avatar
Pedro TOME committed
73 74 75 76 77 78
    # Calculate curvatures
    k = numpy.zeros((img_h, img_w, 4))
    k[:,:,0] = (fxx/((1 + fx**2)**(3/2)))*finger_mask  # hor #
    k[:,:,1] = (fyy/((1 + fy**2)**(3/2)))*finger_mask  # ver #
    k[:,:,2] = (f11/((1 + f1**2)**(3/2)))*finger_mask  # \   #
    k[:,:,3] = (f22/((1 + f2**2)**(3/2)))*finger_mask  # /   #
79

Pedro TOME's avatar
Pedro TOME committed
80 81 82
    # Scores
    Vt = numpy.zeros(image.shape)
    Wr = 0
83

Pedro TOME's avatar
Pedro TOME committed
84 85
    # Horizontal direction
    bla = k[:,:,0] > 0
86 87
    for y in range(0,img_h):
        for x in range(0,img_w):
Pedro TOME's avatar
Pedro TOME committed
88 89 90 91 92 93 94 95
            if (bla[y,x]):
                Wr = Wr + 1
            if ( Wr > 0 and (x == (img_w-1) or not bla[y,x]) ):
                if (x == (img_w-1)):
                    # Reached edge of image
                    pos_end = x
                else:
                    pos_end = x - 1
96 97

                pos_start = pos_end - Wr + 1 # Start pos of concave
Pedro TOME's avatar
Pedro TOME committed
98 99
                if (pos_start == pos_end):
                    I=numpy.argmax(k[y,pos_start,0])
100
                else:
Pedro TOME's avatar
Pedro TOME committed
101
                    I=numpy.argmax(k[y,pos_start:pos_end+1,0])
102

Pedro TOME's avatar
Pedro TOME committed
103 104 105
                pos_max = pos_start + I
                Scr = k[y,pos_max,0]*Wr
                Vt[y,pos_max] = Vt[y,pos_max] + Scr
106 107
                Wr = 0

Pedro TOME's avatar
Pedro TOME committed
108 109 110

    # Vertical direction
    bla = k[:,:,1] > 0
111 112
    for x in range(0,img_w):
        for y in range(0,img_h):
Pedro TOME's avatar
Pedro TOME committed
113 114 115 116 117 118 119
            if (bla[y,x]):
                Wr = Wr + 1
            if ( Wr > 0 and (y == (img_h-1) or not bla[y,x]) ):
                if (y == (img_h-1)):
                    # Reached edge of image
                    pos_end = y
                else:
120 121
                    pos_end = y - 1

Pedro TOME's avatar
Pedro TOME committed
122 123 124 125 126
                pos_start = pos_end - Wr + 1 # Start pos of concave
                if (pos_start == pos_end):
                    I=numpy.argmax(k[pos_start,x,1])
                else:
                    I=numpy.argmax(k[pos_start:pos_end+1,x,1])
127 128

                pos_max = pos_start + I
Pedro TOME's avatar
Pedro TOME committed
129
                Scr = k[pos_max,x,1]*Wr
130

Pedro TOME's avatar
Pedro TOME committed
131 132
                Vt[pos_max,x] = Vt[pos_max,x] + Scr
                Wr = 0
133

Pedro TOME's avatar
Pedro TOME committed
134 135 136 137 138 139 140 141 142 143 144
    # Direction: \ #
    bla = k[:,:,2] > 0
    for start in range(0,img_w+img_h-1):
        # Initial values
        if (start <= img_w-1):
            x = start
            y = 0
        else:
            x = 0
            y = start - img_w + 1
        done = False
145

Pedro TOME's avatar
Pedro TOME committed
146 147 148
        while (not done):
            if(bla[y,x]):
                Wr = Wr + 1
149

Pedro TOME's avatar
Pedro TOME committed
150 151 152 153 154 155 156 157
            if ( Wr > 0 and (y == img_h-1 or x == img_w-1 or not bla[y,x]) ):
                if (y == img_h-1 or x == img_w-1):
                    # Reached edge of image
                    pos_x_end = x
                    pos_y_end = y
                else:
                    pos_x_end = x - 1
                    pos_y_end = y - 1
158

Pedro TOME's avatar
Pedro TOME committed
159 160
                pos_x_start = pos_x_end - Wr + 1
                pos_y_start = pos_y_end - Wr + 1
161

Pedro TOME's avatar
Pedro TOME committed
162 163 164 165 166 167 168 169
                if (pos_y_start == pos_y_end and pos_x_start == pos_x_end):
                    d = k[pos_y_start, pos_x_start, 2]
                elif (pos_y_start == pos_y_end):
                    d = numpy.diag(k[pos_y_start, pos_x_start:pos_x_end+1, 2])
                elif (pos_x_start == pos_x_end):
                    d = numpy.diag(k[pos_y_start:pos_y_end+1, pos_x_start, 2])
                else:
                    d = numpy.diag(k[pos_y_start:pos_y_end+1, pos_x_start:pos_x_end+1, 2])
170 171 172 173 174 175

                I = numpy.argmax(d)

                pos_x_max = pos_x_start + I
                pos_y_max = pos_y_start + I

Pedro TOME's avatar
Pedro TOME committed
176
                Scr = k[pos_y_max,pos_x_max,2]*Wr
177

Pedro TOME's avatar
Pedro TOME committed
178 179
                Vt[pos_y_max,pos_x_max] = Vt[pos_y_max,pos_x_max] + Scr
                Wr = 0
180

Pedro TOME's avatar
Pedro TOME committed
181 182 183 184 185
            if((x == img_w-1) or (y == img_h-1)):
                done = True
            else:
                x = x + 1
                y = y + 1
186

Pedro TOME's avatar
Pedro TOME committed
187 188 189 190 191 192 193 194 195 196 197
    # Direction: /
    bla = k[:,:,3] > 0
    for start in range(0,img_w+img_h-1):
        # Initial values
        if (start <= (img_w-1)):
            x = start
            y = img_h-1
        else:
            x = 0
            y = img_w+img_h-start-1
        done = False
198

Pedro TOME's avatar
Pedro TOME committed
199 200 201 202 203 204 205 206 207 208 209
        while (not done):
            if(bla[y,x]):
                Wr = Wr + 1
            if ( Wr > 0 and (y == 0 or x == img_w-1 or not bla[y,x]) ):
                if (y == 0 or x == img_w-1):
                    # Reached edge of image
                    pos_x_end = x
                    pos_y_end = y
                else:
                    pos_x_end = x - 1
                    pos_y_end = y + 1
210

Pedro TOME's avatar
Pedro TOME committed
211 212
                pos_x_start = pos_x_end - Wr + 1
                pos_y_start = pos_y_end + Wr - 1
213

Pedro TOME's avatar
Pedro TOME committed
214 215 216 217 218 219 220 221
                if (pos_y_start == pos_y_end and pos_x_start == pos_x_end):
                    d = k[pos_y_end, pos_x_start, 3]
                elif (pos_y_start == pos_y_end):
                    d = numpy.diag(numpy.flipud(k[pos_y_end, pos_x_start:pos_x_end+1, 3]))
                elif (pos_x_start == pos_x_end):
                    d = numpy.diag(numpy.flipud(k[pos_y_end:pos_y_start+1, pos_x_start, 3]))
                else:
                    d = numpy.diag(numpy.flipud(k[pos_y_end:pos_y_start+1, pos_x_start:pos_x_end+1, 3]))
222 223 224 225

                I = numpy.argmax(d)
                pos_x_max = pos_x_start + I
                pos_y_max = pos_y_start - I
Pedro TOME's avatar
Pedro TOME committed
226 227 228
                Scr = k[pos_y_max,pos_x_max,3]*Wr
                Vt[pos_y_max,pos_x_max] = Vt[pos_y_max,pos_x_max] + Scr
                Wr = 0
229

Pedro TOME's avatar
Pedro TOME committed
230 231 232 233
            if((x == img_w-1) or (y == 0)):
                done = True
            else:
                x = x + 1
234
                y = y - 1
Pedro TOME's avatar
Pedro TOME committed
235 236 237 238 239 240 241 242 243

    ## Connection of vein centres
    Cd = numpy.zeros((img_h, img_w, 4))
    for x in range(2,img_w-3):
        for y in range(2,img_h-3):
            Cd[y,x,0] = min(numpy.amax(Vt[y,x+1:x+3]), numpy.amax(Vt[y,x-2:x]))   # Hor  #
            Cd[y,x,1] = min(numpy.amax(Vt[y+1:y+3,x]), numpy.amax(Vt[y-2:y,x]))   # Vert #
            Cd[y,x,2] = min(numpy.amax(Vt[y-2:y,x-2:x]), numpy.amax(Vt[y+1:y+3,x+1:x+3])) # \  #
            Cd[y,x,3] = min(numpy.amax(Vt[y+1:y+3,x-2:x]), numpy.amax(Vt[y-2:y,x+1:x+3])) # /  #
244

Pedro TOME's avatar
Pedro TOME committed
245 246
    #Veins
    img_veins = numpy.amax(Cd,axis=2)
247

Pedro TOME's avatar
Pedro TOME committed
248 249 250 251 252
    # Binarise the vein image
    md = numpy.median(img_veins[img_veins>0])
    img_veins_bin = img_veins > md

    return img_veins_bin.astype(numpy.float64)
253 254 255


  def __call__(self, image):
Pedro TOME's avatar
Pedro TOME committed
256
    """Reads the input image, extract the features based on Maximum Curvature of the fingervein image, and writes the resulting template"""
257

258
    finger_image = image[0] #Normalized image with or without histogram equalization
259 260 261
    finger_mask = image[1]

    return self.maximum_curvature(finger_image, finger_mask)