MaximumCurvature.py 8.78 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
11
from bob.bio.base.features.Extractor import Extractor

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
  """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
Pedro TOME's avatar
Pedro TOME committed
21
  """
22
23


Pedro TOME's avatar
Pedro TOME committed
24
25
26
27
28
29
30
31
32
33
34
35
  def __init__(
      self,
      sigma = 5, #Sigma used for determining derivatives
      gpu = False
  ):

    # call base class constructor
    Extractor.__init__(
        self,
        sigma = sigma,
        gpu = gpu
    )
36

Pedro TOME's avatar
Pedro TOME committed
37
38
39
    # block parameters
    self.sigma = sigma
    self.gpu = gpu
40
41
42
43
44
45


  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
46
47
48
49
50
    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)
51
52
    finger_mask[mask == True] = 1

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

Pedro TOME's avatar
Pedro TOME committed
55
56
57
    x = numpy.arange(-winsize, winsize+1)
    y = numpy.arange(-winsize, winsize+1)
    X, Y = numpy.meshgrid(x, y)
58

Pedro TOME's avatar
Pedro TOME committed
59
60
61
62
63
64
    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
65

Pedro TOME's avatar
Pedro TOME committed
66
    # Do the actual filtering
67

Pedro TOME's avatar
Pedro TOME committed
68
69
70
71
72
    fx = utils.imfilter(image, hx, self.gpu, conv=False)
    fxx = utils.imfilter(image, hxx, self.gpu, conv=False)
    fy = utils.imfilter(image, hy, self.gpu, conv=False)
    fyy = utils.imfilter(image, hyy, self.gpu, conv=False)
    fxy = utils.imfilter(image, hxy, self.gpu, conv=False)
73

Pedro TOME's avatar
Pedro TOME committed
74
75
76
77
    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       # // #
78

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

Pedro TOME's avatar
Pedro TOME committed
81
82
83
84
85
86
    # 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  # /   #
87

Pedro TOME's avatar
Pedro TOME committed
88
89
90
    # Scores
    Vt = numpy.zeros(image.shape)
    Wr = 0
91

Pedro TOME's avatar
Pedro TOME committed
92
93
    # Horizontal direction
    bla = k[:,:,0] > 0
94
95
    for y in range(0,img_h):
        for x in range(0,img_w):
Pedro TOME's avatar
Pedro TOME committed
96
97
98
99
100
101
102
103
            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
104
105

                pos_start = pos_end - Wr + 1 # Start pos of concave
Pedro TOME's avatar
Pedro TOME committed
106
107
                if (pos_start == pos_end):
                    I=numpy.argmax(k[y,pos_start,0])
108
                else:
Pedro TOME's avatar
Pedro TOME committed
109
                    I=numpy.argmax(k[y,pos_start:pos_end+1,0])
110

Pedro TOME's avatar
Pedro TOME committed
111
112
113
                pos_max = pos_start + I
                Scr = k[y,pos_max,0]*Wr
                Vt[y,pos_max] = Vt[y,pos_max] + Scr
114
115
                Wr = 0

Pedro TOME's avatar
Pedro TOME committed
116
117
118

    # Vertical direction
    bla = k[:,:,1] > 0
119
120
    for x in range(0,img_w):
        for y in range(0,img_h):
Pedro TOME's avatar
Pedro TOME committed
121
122
123
124
125
126
127
            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:
128
129
                    pos_end = y - 1

Pedro TOME's avatar
Pedro TOME committed
130
131
132
133
134
                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])
135
136

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

Pedro TOME's avatar
Pedro TOME committed
139
140
                Vt[pos_max,x] = Vt[pos_max,x] + Scr
                Wr = 0
141

Pedro TOME's avatar
Pedro TOME committed
142
143
144
145
146
147
148
149
150
151
152
    # 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
153

Pedro TOME's avatar
Pedro TOME committed
154
155
156
        while (not done):
            if(bla[y,x]):
                Wr = Wr + 1
157

Pedro TOME's avatar
Pedro TOME committed
158
159
160
161
162
163
164
165
            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
166

Pedro TOME's avatar
Pedro TOME committed
167
168
                pos_x_start = pos_x_end - Wr + 1
                pos_y_start = pos_y_end - Wr + 1
169

Pedro TOME's avatar
Pedro TOME committed
170
171
172
173
174
175
176
177
                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])
178
179
180
181
182
183

                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
184
                Scr = k[pos_y_max,pos_x_max,2]*Wr
185

Pedro TOME's avatar
Pedro TOME committed
186
187
                Vt[pos_y_max,pos_x_max] = Vt[pos_y_max,pos_x_max] + Scr
                Wr = 0
188

Pedro TOME's avatar
Pedro TOME committed
189
190
191
192
193
            if((x == img_w-1) or (y == img_h-1)):
                done = True
            else:
                x = x + 1
                y = y + 1
194

Pedro TOME's avatar
Pedro TOME committed
195
196
197
198
199
200
201
202
203
204
205
    # 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
206

Pedro TOME's avatar
Pedro TOME committed
207
208
209
210
211
212
213
214
215
216
217
        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
218

Pedro TOME's avatar
Pedro TOME committed
219
220
                pos_x_start = pos_x_end - Wr + 1
                pos_y_start = pos_y_end + Wr - 1
221

Pedro TOME's avatar
Pedro TOME committed
222
223
224
225
226
227
228
229
                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]))
230
231
232
233

                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
234
235
236
                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
237

Pedro TOME's avatar
Pedro TOME committed
238
239
240
241
            if((x == img_w-1) or (y == 0)):
                done = True
            else:
                x = x + 1
242
                y = y - 1
Pedro TOME's avatar
Pedro TOME committed
243
244
245
246
247
248
249
250
251

    ## 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])) # /  #
252

Pedro TOME's avatar
Pedro TOME committed
253
254
    #Veins
    img_veins = numpy.amax(Cd,axis=2)
255

Pedro TOME's avatar
Pedro TOME committed
256
257
258
259
260
    # 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)
261
262
263


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

Pedro TOME's avatar
Pedro TOME committed
266
    finger_image = image[0]    #Normalized image with or without histogram equalization
267
268
269
270
    finger_mask = image[1]

    return self.maximum_curvature(finger_image, finger_mask)

Pedro TOME's avatar
Pedro TOME committed
271
272
273
274

  def save_feature(self, feature, feature_file):
    f = bob.io.base.HDF5File(feature_file, 'w')
    f.set('feature', feature)
275

Pedro TOME's avatar
Pedro TOME committed
276
277
278
  def read_feature(self, feature_file):
    f = bob.io.base.HDF5File(feature_file, 'r')
    image = f.read('feature')
279
    return image