Dataset:https://personalpages.manchester.ac.uk/staff/timothy.f.cootes/data/talking_face/talking_face.html
step 1: load images and points
import cv2
from matplotlib import pyplot as plt
import os
import math
import numpy as np
%matplotlib inline
random_ind = np.random.randint(0, 4999, size=(50,))
images_dir = "./Talking Face Video/images"
images = []
for i in random_ind:
file_name = "franck_"+"{:05d}".format(i)+".jpg"
images.append(cv2.cvtColor(cv2.imread(os.path.join(images_dir, file_name)), cv2.COLOR_RGB2BGR))
points_dir = "./Talking Face Video/points"
points_images = []
for i in random_ind:
file_name = "franck_"+"{:05d}".format(i)+".pts"
with open(os.path.join(points_dir, file_name),'r') as f:
lines = f.readlines()
points_images.append( [tuple([int(float(num)) for num in line.split()]) for line in lines[3:71]])
step 2: select the ROI of faces
faceCascade = cv2.CascadeClassifier("/usr/local/Cellar/opencv/4.0.1/share/opencv4/haarcascades/haarcascade_frontalface_default.xml")
faces = []
points_faces = []
for image, points_image in zip(images, points_images):
fcs = faceCascade.detectMultiScale(
image,
scaleFactor=1.1,
minNeighbors=5,
minSize=(30, 30),
flags = 0
)
for (x, y, w, h) in fcs:
fcs_img = image[y:y+h+70, x:x+w]
points_face = [tuple(np.subtract(pt, (x,y))) for pt in points_image]
faces.append(fcs_img)
points_faces.append(points_face)
plot function defininitions
def get_landmarks_image(image, points, labeled=True, point_size=2):
image_landmarks = image.copy()
ind = 0
for point in points:
cv2.circle(image_landmarks, point, point_size, (255,0,0), -1)
if labeled:
cv2.putText(image_landmarks, str(ind), tuple(np.add(point, (5,0))), cv2.FONT_HERSHEY_SIMPLEX, 0.3, (0,255,0), 1, cv2.LINE_AA)
ind = ind+1
return image_landmarks
def get_connected_landmarks_image(image, points, line_color=(0,0,255), point_size=2, line_size=1):
image_landmarks = image.copy()
for ind_points in range(len(points)):
# plot the contour of face
if 0 < ind_points < 15:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
# plot the contour of right eyebrow
elif 15 < ind_points < 21:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 21:
cv2.line(image_landmarks, points[20], points[15], line_color, line_size)
# plot the contour of left eyebrow
elif 21 < ind_points < 27:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 27:
cv2.line(image_landmarks, points[26], points[21], line_color, line_size)
# plot the contour of left eye
elif 27 < ind_points < 31:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 31:
cv2.line(image_landmarks, points[30], points[27], line_color, line_size)
# plot the contour of right eye
elif 32 < ind_points < 36:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 36:
cv2.line(image_landmarks, points[35], points[32], line_color, line_size)
# plot the contour of nose
elif 37 < ind_points < 46:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 46 or ind_points == 47:
cv2.line(image_landmarks, points[ind_points], points[ind_points-6], line_color, line_size)
cv2.line(image_landmarks, points[ind_points], points[ind_points-5], line_color, line_size)
# plot the contour of external lips
elif 48 < ind_points < 60:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 60:
cv2.line(image_landmarks, points[59], points[48], line_color, line_size)
cv2.line(image_landmarks, points[60], points[48], line_color, line_size)
# plot the contour of internal lips
elif 60 < ind_points < 63:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 63:
cv2.line(image_landmarks, points[62], points[54], line_color, line_size)
cv2.line(image_landmarks, points[63], points[54], line_color, line_size)
elif 63 < ind_points < 66:
cv2.line(image_landmarks, points[ind_points], points[ind_points-1], line_color, line_size)
elif ind_points == 66:
cv2.line(image_landmarks, points[65], points[48], line_color, line_size)
cv2.circle(image_landmarks, points[ind_points], point_size, (255,0,0), -1)
return image_landmarks
def get_points_on_perpendicular_line(intersection_point, pt1, pt2, length):
rotation_matrix = np.array([[math.cos(math.pi/2), -math.sin(math.pi/2)],
[math.sin(math.pi/2), math.cos(math.pi/2)]])
normal = np.dot(rotation_matrix, np.subtract(pt2, pt1))
pt1_normal = np.add(intersection_point, normal/np.sqrt(normal.dot(normal))*(length/2)).astype('int64')
pt2_normal = np.subtract(intersection_point, normal/np.sqrt(normal.dot(normal))*(length/2)).astype('int64')
return [tuple(pt1_normal), tuple(pt2_normal)]
def get_landmarks_normal_image(image, points):
image_landmarks = image.copy()
for ind_points in range(len(points)):
if ind_points != 14 and ind_points !=31 and ind_points !=36 and ind_points !=45 and ind_points != 66 and ind_points != 67:
normal_points = get_points_on_perpendicular_line(points[ind_points], points[ind_points], points[ind_points+1], 20)
cv2.line(image_landmarks, normal_points[0], points[ind_points], (0,0,255), 1)
cv2.line(image_landmarks, points[ind_points], normal_points[1], (0,255,0), 1)
elif ind_points == 14 or ind_points == 45:
normal_points = get_points_on_perpendicular_line(points[ind_points], points[ind_points-1], points[ind_points], 20)
cv2.line(image_landmarks, normal_points[0], points[ind_points], (0,0,255), 1)
cv2.line(image_landmarks, points[ind_points], normal_points[1], (0,255,0), 1)
cv2.circle(image_landmarks, points[ind_points], 2, (255,0,0), -1)
return image_landmarks
plt.imshow(get_landmarks_normal_image(faces[0], points_faces[0]))
<matplotlib.image.AxesImage at 0x11ee2ec18>
procustes transform
def get_mean_shape_landmarks(points):
mean_shape = points[0].copy()
for i in range(1,len(points)):
for j in range(len(mean_shape)):
mean_shape[j] = tuple(np.add(mean_shape[j], points[i][j]))
return [tuple(np.divide(i, len(points)).astype('int64')) for i in mean_shape]
def get_centroid_size(points):
centroid_point = (np.sum(np.array(points), axis=0)/len(points)).astype('int')
vectors_from_centroid = np.array(points) - centroid_point
centroid_size = 0
for vector in vectors_from_centroid:
centroid_size = centroid_size + vector.dot(vector)
return np.sqrt(centroid_size)
def scaling(pts_reference, pts):
centroid_size_pts_reference = get_centroid_size(pts_reference)
centroid_size_pts = get_centroid_size(pts)
ratio = centroid_size_pts_reference/centroid_size_pts
return [tuple(i) for i in np.multiply(pts, ratio).astype('int')], ratio
def translating(pts_reference, pts):
centroid_point_pts_reference = (np.sum(np.array(pts_reference), axis=0)/len(pts_reference)).astype('int')
centroid_point_pts = (np.sum(np.array(pts), axis=0)/len(pts)).astype('int')
translation = np.subtract(centroid_point_pts_reference, centroid_point_pts)
return [tuple(i) for i in np.add(np.array(pts), translation)], translation
def rotating(pts_reference, pts, degree):
degree = degree*math.pi/180
rotation_matrix = np.array([[math.cos(degree), -math.sin(degree)],
[math.sin(degree), math.cos(degree)]])
pts_rotated, _ = translating(pts_reference, [tuple(np.dot(rotation_matrix, pt).astype('int')) for pt in pts])
return pts_rotated
def get_euclidean_distance_between_shapes(pts1, pts2):
euclidean_distance = 0
for pt1,pt2 in zip(pts1, pts2):
difference_vector = np.subtract(pt1, pt2)
euclidean_distance = euclidean_distance + difference_vector.dot(difference_vector)
return np.sqrt(euclidean_distance)
def procustes_transform(pts_reference, pts):
pts_scaled, ratio = scaling(pts_reference, pts)
pts_scaled_translated, translation = translating(pts_reference, pts_scaled)
degree_minDist = 360
minDist = np.finfo('float').max
for degree in np.linspace(-15*math.pi/180, 15*math.pi/180, 100):
dist = get_euclidean_distance_between_shapes(pts_reference, rotating(pts_reference, pts_scaled_translated, degree))
if dist < minDist:
minDist = dist
degree_minDist = degree
return rotating(pts_reference, pts_scaled_translated, degree_minDist), [ratio, degree_minDist, translation]
def get_local_area(face, points, radius):
local_areas = []
for pt in points:
local_areas.append(face[pt[1]-radius:pt[1]+radius, pt[0]-radius:pt[0]+radius])
return local_areas
def get_local_areas_faces(faces, points_faces, radius):
local_areas_faces = []
for face,points_face in zip(faces,points_faces):
local_areas_faces.append(get_local_area(face, points_face, radius))
return local_areas_faces
def get_local_areas_points(local_areas_faces):
local_areas_points = []
for i in range(len(local_areas_faces[0])):
local_areas_point = []
for j in range(len(local_areas_faces)):
local_areas_point.append(local_areas_faces[j][i])
local_areas_points.append(local_areas_point)
return local_areas_points
def centering(points, centroid):
translating_vector = np.subtract(centroid, points[67])
return [tuple(i) for i in np.add(np.array(points), translating_vector)]
def scaling_ratio(points, ratio):
return [tuple(i) for i in np.multiply(points, ratio).astype('int')]
PCA analysis
def get_principal_components(trainning_X, num_of_components):
mean_x = np.mean(trainning_X, axis=1).reshape(-1,1)
normalized_X = trainning_X - mean_x
covariance_matrix = np.dot(normalized_X.astype('float').transpose(), normalized_X.astype('float'))/normalized_X.shape[1]
eigenvalues, eigenvectors = np.linalg.eig(covariance_matrix)
eigenvectors = np.dot(normalized_X.astype('float'), eigenvectors)
principal_components = eigenvectors[:,np.argsort(eigenvalues)[-num_of_components]].reshape(-1,1)
for i in np.argsort(eigenvalues)[-(num_of_components-1):]:
principal_components = np.concatenate((principal_components, eigenvectors[:,i].reshape(-1,1)), axis=1)
''''
eigenfaces_matrix = 255*(principal_components-np.min(principal_components))/(np.max(principal_components)-np.min(principal_components))
eigenfaces = []
for i in range(eigenfaces_matrix.shape[1]):
eigenfaces.append(eigenfaces_matrix[:,i].astype('uint8').reshape(images[0].shape[0], images[0].shape[1]))
eigenfaces.reverse()
plt.figure(figsize=(20,10))
for i in range(len(eigenfaces)):
plt.subplot(5,10,i+1)
plt.imshow(eigenfaces[i], cmap=plt.cm.gray)
'''
for i in range(principal_components.shape[1]):
eigenvector_length = math.sqrt(np.dot(principal_components[:,i].transpose(), principal_components[:,i]))
principal_components[:,i] = principal_components[:,i]/eigenvector_length
return principal_components, mean_x, np.sort(eigenvalues)[-num_of_components:]
def get_min_dist(local_areas_PCA_space, vector_to_compare_PCA_space):
minDist = np.finfo('float').max
for i in range(local_areas_PCA_space.shape[1]):
dist = np.sqrt(np.sum(np.power(local_areas_PCA_space[:,i]-vector_to_compare_PCA_space,2)))
if dist < minDist:
minDist = dist
return minDist
def get_best_position(search_area, local_areas_PCA_space,principal_components, mean_local_areas,
point, minDist_point, radius, search_radius, search_step):
minDist = np.finfo('float').max
best_position = (0,0)
best_area = np.zeros((2*radius,2*radius))
for x in range(0, search_area.shape[1]-radius*2+1, search_step):
for y in range(0, search_area.shape[0]-radius*2+1, search_step):
area_to_compare = search_area[y:y+2*radius, x:x+2*radius]
vector_to_compare = area_to_compare.reshape(-1,1)
vector_to_compare_PCA_space = np.dot(principal_components.transpose(), (vector_to_compare-mean_local_areas).astype('float'))
dist = get_min_dist(local_areas_PCA_space, vector_to_compare_PCA_space)
if dist < minDist:
minDist = dist
best_position = (x+radius,y+radius)
best_area = area_to_compare
if minDist < minDist_point:
return minDist, tuple(np.add(np.subtract(point, (search_radius, search_radius)), best_position)), best_area, best_position
else:
return minDist_point, point, best_area, best_position
def ASM_search(face_gray, local_areas_face_gray, asm_model, points_PCA_space, points_principal_components, points_mean_areas,
minDist_asm_model, radius, search_radius, search_step):
minDist_asm_model = np.ones((len(asm_model), 1), dtype='float')*np.finfo('float').max
tmp_asm_model = asm_model.copy()
for j in range(len(tmp_asm_model)):
# for j in range(27,32):
search_area = face_gray[tmp_asm_model[j][1]-search_radius:tmp_asm_model[j][1]+search_radius,
tmp_asm_model[j][0]-search_radius:tmp_asm_model[j][0]+search_radius]
# plt.figure(figsize=(8,3))
# plt.suptitle("landmark {}: {}".format(j,tmp_asm_model[j]))
minDist_asm_model[j],tmp_asm_model[j],best_area,best_position = get_best_position(search_area, points_PCA_space[j], points_principal_components[j],
points_mean_areas[j], tmp_asm_model[j], minDist_asm_model[j], radius, search_radius, search_step)
# plt.subplot(1,3,1)
# plt.imshow(search_area, cmap=plt.cm.gray)
# plt.title("-{}+{}={}".format((search_radius,search_radius), best_position, tmp_asm_model[j]))
# plt.subplot(1,3,2)
# plt.imshow(local_areas_face_gray[j], cmap=plt.cm.gray)
# plt.title("target area")
# plt.subplot(1,3,3)
# plt.imshow(best_area, cmap=plt.cm.gray)
# plt.title(str(minDist_asm_model[j]))
return tmp_asm_model, minDist_asm_model
image pyramid level 1
centroid_point = (350,350)
figshape = (700,700,3)
search_step = 1
iteration = 5
faces_L1 = [cv2.pyrDown(face) for face in faces]
points_faces_L1 = []
for points_face in points_faces:
points_face_L1 = []
for point in points_face:
points_face_L1.append((int(point[0]/2), int(point[1]/2)))
points_faces_L1.append(points_face_L1)
radius = np.min(points_faces_L1)
search_radius = 3*radius
mean_shape_L1 = get_mean_shape_landmarks(points_faces_L1)
points_faces_aligned_L1 = [procustes_transform(mean_shape_L1, pts)[0] for pts in points_faces_L1]
mean_shape_aligned_L1 = get_mean_shape_landmarks(points_faces_aligned_L1)
asm_model_L1 = centering(scaling_ratio(mean_shape_aligned_L1, 0.8), (int(faces_L1[0].shape[1]/2), int(faces_L1[0].shape[0]/2)))
asm_model_L1 = [(int(i[0]), int(i[1])) for i in asm_model_L1]
faces_gray_L1 = [cv2.cvtColor(face, cv2.COLOR_BGR2GRAY) for face in faces_L1]
local_areas_faces_gray_L1 = get_local_areas_faces(faces_gray_L1, points_faces_L1, radius)
local_areas_gray_points_L1 = get_local_areas_points(local_areas_faces_gray_L1)
trainning_shapes_L1 = np.array(points_faces_L1[0]).reshape(-1,1)
for i in range(1,len(points_faces_L1)):
trainning_shapes_L1 = np.hstack((trainning_shapes_L1, np.array(points_faces_L1[i]).reshape(-1,1)))
asm_model_shape_L1 = np.array(asm_model_L1).reshape(-1,1)
shape_principal_components_L1, mean_trainning_shape_L1, eigenvalues_L1 = get_principal_components(trainning_shapes_L1, 8)
local_areas_points_X_L1 = []
for i in range(len(local_areas_gray_points_L1)):
local_areas_point_X_L1 = np.array(local_areas_gray_points_L1[i][0]).reshape(-1,1)
for j in range(1,len(local_areas_gray_points_L1[i])):
local_areas_point_X_L1 = np.hstack((local_areas_point_X_L1, np.array(local_areas_gray_points_L1[i][j]).reshape(-1,1)))
local_areas_points_X_L1.append(local_areas_point_X_L1)
local_areas_points_principal_components_L1 = [get_principal_components(local_areas_point_X_L1, 8)[0]
for local_areas_point_X_L1 in local_areas_points_X_L1]
mean_local_areas_points_L1 = [get_principal_components(local_areas_point_X_L1, 8)[1] for local_areas_point_X_L1 in local_areas_points_X_L1]
local_areas_points_PCA_space_L1 = [np.dot(i.transpose(), (j-h).astype('float')) for i,j,h in
zip(local_areas_points_principal_components_L1, local_areas_points_X_L1, mean_local_areas_points_L1)]
asm_model_PCA_space_L1 = np.dot(shape_principal_components_L1.transpose(), np.subtract(asm_model_shape_L1, mean_trainning_shape_L1))
PCA_space_range = np.add([(-3*np.sqrt(eigenvalue), 3*np.sqrt(eigenvalue)) for eigenvalue in eigenvalues_L1], asm_model_PCA_space_L1)
ASM search showed by connected landmarks
plt.figure(figsize=(20,5))
plt.subplot(1,6,1)
plt.imshow(get_connected_landmarks_image(faces_L1[0], asm_model_L1, point_size=1, line_size=1))
plt.title("initial ASM")
minDist_asm_model_L1 = np.ones((len(asm_model_L1), 1), dtype='float')*np.finfo('float').max
for i in range(iteration):
asm_model_L1_new, minDist_asm_model_L1 = ASM_search(faces_gray_L1[0], local_areas_faces_gray_L1[0], asm_model_L1,
local_areas_points_PCA_space_L1, local_areas_points_principal_components_L1,
mean_local_areas_points_L1, minDist_asm_model_L1, radius, search_radius, search_step)
asm_model_L1_new_aligned, _ = procustes_transform(asm_model_L1, asm_model_L1_new)
asm_model_L1_new_aligned_PCA_space = np.dot(shape_principal_components_L1.transpose(),
np.subtract(np.array(asm_model_L1_new_aligned).reshape(-1,1), mean_trainning_shape_L1))
for j in range(len(asm_model_L1_new_aligned_PCA_space)):
if asm_model_L1_new_aligned_PCA_space[j][0] > PCA_space_range[j][1]:
asm_model_L1_new_aligned_PCA_space[j][0] = PCA_space_range[j][1]
elif asm_model_L1_new_aligned_PCA_space[j][0] < PCA_space_range[j][0]:
asm_model_L1_new_aligned_PCA_space[j][0] = PCA_space_range[j][0]
asm_model_L1_new_aligned_recovered = np.dot(shape_principal_components_L1, asm_model_L1_new_aligned_PCA_space) + mean_trainning_shape_L1
asm_model_L1_new_aligned_recovered = [(int(asm_model_L1_new_aligned_recovered[i][0]), int(asm_model_L1_new_aligned_recovered[i+1][0]))
for i in range(0, asm_model_L1_new_aligned_recovered.shape[0], 2)]
asm_model_L1, _ = procustes_transform(asm_model_L1_new, asm_model_L1_new_aligned_recovered)
plt.subplot(1,6,i+2)
plt.imshow(get_connected_landmarks_image(faces_L1[0], asm_model_L1, point_size=1, line_size=1))
plt.title("ASM after {} times iteration".format(i+1))
ASM search showed by triangles patches
def is_inside(rect, point):
if rect[0] <= point[0] <= rect[0]+rect[2] and rect[1] <= point[1] <= rect[1]+rect[3]:
return True
else:
return False
def triangulating(image, points, line_size=1, point_size=3):
face = image.copy()
convexhull = cv2.convexHull(np.array(points))
points_inds = {points[i]:i for i in range(len(points))}
rect = cv2.boundingRect(convexhull)
subdiv = cv2.Subdiv2D(rect)
subdiv.insert(points)
triangles = subdiv.getTriangleList().astype('int')
triangles_inds = []
for triangle in triangles:
pt1 = (triangle[0], triangle[1])
pt2 = (triangle[2], triangle[3])
pt3 = (triangle[4], triangle[5])
if is_inside(rect, pt1) and is_inside(rect, pt2) and is_inside(rect, pt3):
cv2.line(face, pt1, pt2, (0,0,255),line_size)
cv2.line(face, pt2, pt3, (0,0,255),line_size)
cv2.line(face, pt3, pt1, (0,0,255),line_size)
triangles_inds.append([points_inds[pt1], points_inds[pt2], points_inds[pt3]])
face = get_landmarks_image(face, points, labeled=False, point_size=point_size)
return face, triangles_inds
plt.figure(figsize=(20,5))
plt.subplot(1,6,1)
# plt.figure(figsize=(3,4))
initial_asm, triangles_inds = triangulating(faces_L1[0], asm_model_L1, 1, 1)
plt.imshow(initial_asm)
plt.title("initial ASM")
# plt.savefig("initial_asm.png")
minDist_asm_model_L1 = np.ones((len(asm_model_L1), 1), dtype='float')*np.finfo('float').max
for i in range(iteration):
asm_model_L1_new, minDist_asm_model_L1 = ASM_search(faces_gray_L1[0], local_areas_faces_gray_L1[0], asm_model_L1,
local_areas_points_PCA_space_L1, local_areas_points_principal_components_L1,
mean_local_areas_points_L1, minDist_asm_model_L1, radius, search_radius, search_step)
asm_model_L1_new_aligned, _ = procustes_transform(asm_model_L1, asm_model_L1_new)
asm_model_L1_new_aligned_PCA_space = np.dot(shape_principal_components_L1.transpose(),
np.subtract(np.array(asm_model_L1_new_aligned).reshape(-1,1), mean_trainning_shape_L1))
for j in range(len(asm_model_L1_new_aligned_PCA_space)):
if asm_model_L1_new_aligned_PCA_space[j][0] > PCA_space_range[j][1]:
asm_model_L1_new_aligned_PCA_space[j][0] = PCA_space_range[j][1]
elif asm_model_L1_new_aligned_PCA_space[j][0] < PCA_space_range[j][0]:
asm_model_L1_new_aligned_PCA_space[j][0] = PCA_space_range[j][0]
asm_model_L1_new_aligned_recovered = np.dot(shape_principal_components_L1, asm_model_L1_new_aligned_PCA_space) + mean_trainning_shape_L1
asm_model_L1_new_aligned_recovered = [(int(asm_model_L1_new_aligned_recovered[i][0]), int(asm_model_L1_new_aligned_recovered[i+1][0]))
for i in range(0, asm_model_L1_new_aligned_recovered.shape[0], 2)]
asm_model_L1, _ = procustes_transform(asm_model_L1_new, asm_model_L1_new_aligned_recovered)
inds_points = {i:asm_model_L1[i] for i in range(len(asm_model_L1))}
convexhull = cv2.convexHull(np.array(asm_model_L1))
tmp_face = faces_L1[0].copy()
for triangle_inds in triangles_inds:
pt1 = inds_points[triangle_inds[0]]
pt2 = inds_points[triangle_inds[1]]
pt3 = inds_points[triangle_inds[2]]
cv2.line(tmp_face, pt1, pt2, (0,0,255),1)
cv2.line(tmp_face, pt2, pt3, (0,0,255),1)
cv2.line(tmp_face, pt3, pt1, (0,0,255),1)
tmp_face = get_landmarks_image(tmp_face, asm_model_L1, False, 1)
plt.subplot(1,6,i+2)
# plt.figure(figsize=(3,4))
plt.imshow(tmp_face)
plt.title("ASM after {} times iteration".format(i+1))
# plt.savefig("asm_iteration{}.png".format(i+1))
扫一扫
1841
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
