본문 바로가기

OpenCV

homography matrix 다루기 / cv2.findHomography( ), cv2.warpPerspective( )

320x100
320x100

 

 

 제목을 뭐로 해야할지 몰라서 homography matrix 다루기라 두루뭉실하게 적어놨다. 오늘은 아래의 그림과 목차 흐름으로 정리하고자 한다.

 

  1. cv2.xfeatures2d.SIFT_create( ) 함수로 keypoint와 descriptor 추출
  2. cv2.BFMatcher( ) 함수로 keypoint 매칭
  3. cv2.findHomography( ) 함수로 homography matrix, H 추출
  4. 추출한 H로 이미지 변환
    → cv2.warpPerspective( ) 함수 사용하면 쉽게 할 수 있는데,
       이번 글에서는 H를 직접 다루는 방법으로 정리하고자 한다.
  5. 정합 결과 확인

 

 먼저 말해두자면, 난 SIFT 이론에 대해 아직은 정확히 모른다. -> 20.01.15 논문 쓰려고 공부했다. 이해한 데까지 정리도 했다. [openCV] - SIFT 알고리즘 여기를 보면 된다!

 

 homography matrix 이론에 대해서도 정확하게 아는 건 아니지만 개발에 필요해서 구현 방법을 중심으로 짧게 나마 알아봤고, 그 내용들에 대해 정리해놓고자 한다.

 

 SIFT 사용법과 matching 방법은 앞선 글 cv2.BFMatcher( )에 대해 알아보자 에 설명해뒀으니 넘어가도록 하겠다.

 

 


 

 

1. Homography Matrix, H

 3x3 행렬로 변환 행렬에 해당되는 H는 아래와 같이 표현되며 cv2.findHomography( ) 함수를 통해 구해줄 수 있다.

 

 

 

 image A에서 뽑은 keypoint와 매칭되는 image B의 keypoint를 cv2.findHomography( ) 함수에 넣어주면 된다. 이렇게 구한 H를 사용하면 우리는 image A와 image B 정합시켜줄 수 있다. 더 정확히 표현하면 image A를 image B에 겹치게 변환시켜줄 수 있다. 설명이 헷갈릴텐데 밑에 코드와 결과를 같이 보면 충분히 이해할 수 있을 거다.

 

 변환 과정에는 H와 cv2.warpPerspective( ) 함수를 사용하면 되는데 이번 글에선 OpenCV 함수를 안 쓰고 직접 구현하는 방법으로 정리하려고 한다. ( cv2.warpPerspective( ) 함수 사용법은 구글링 하면 많이 나온다. )

 

 

 방법은 생각보다 간단하다. H와 좌표계를 행렬곱 연산을 해주면 된다. 위의 수식은 homography matrix를 검색해보면 자주 검색되는 식인데, 중요한 점은 위의 식 그대로 연산을 진행해주면 안 되고 아래 식으로 계산해야한다는 점이다.

 

 

x'를 정확하게 구하려면

h11*x    +    h12*y    +    h13 에서 계산을 끝낼 게 아니라

이를 h31*x    + h32*y    +    1로 나눠줘야 한다.

그래야 올바른 결과를 얻을 수 있다.

 

이론 설명은 이쯤에서 멈추고 코드를 통해 이해해보도록 하자.

 

 

 


 

 

< 왼쪽 사진 >
< 오른쪽 사진 >
import numpy as np
import cv2


imageA = cv2.imread('./pano_2.jpeg') # 오른쪽 사진
imageB = cv2.imread('./pano_1.jpeg') # 왼쪽 사진

grayA = cv2.cvtColor(imageA,cv2.COLOR_BGR2GRAY)
grayB = cv2.cvtColor(imageB,cv2.COLOR_BGR2GRAY)


sift = cv2.xfeatures2d.SIFT_create()
kpA, desA = sift.detectAndCompute(grayA, None)
kpB, desB = sift.detectAndCompute(grayB, None)


bf = cv2.BFMatcher()
matches = bf.match(desA, desB)


sorted_matches = sorted(matches, key = lambda x : x.distance)
res = cv2.drawMatches(imageA, kpA, imageB, kpB, sorted_matches[:30], None, flags = 2)


src = np.float32([kpA[m.queryIdx].pt for m in matches]).reshape((-1, 1, 2))
dst = np.float32([kpB[m.trainIdx].pt for m in matches]).reshape((-1, 1, 2))
H, status = cv2.findHomography(src, dst, cv2.RANSAC, 5.0)


before = []
for x in range(imageA.shape[1]):
    for y in range(imageA.shape[0]):
        point = [x, y, 1]
        before.append(point)
before = np.array(before).transpose()

after = np.matmul(H, before)
after = after / after[2, :]
after = after[:2, :]
after = np.round(after, 0).astype(np.int)


height, width, _ = imageA.shape
result = np.zeros((height, width * 2, 3), dtype = np.uint8)
for pt1, pt2 in zip(before[:2, :].transpose(), after.transpose()):
    if pt2[1] >= height:
        continue

    if np.sum(pt2 < 0) >= 1:
        continue
    
    result[pt2[1], pt2[0]] = imageA[pt1[1], pt1[0]]
result[0: height, 0 : width] = imageB

cv2.imshow('result', result)
  • [Line 5, 6]
    우리는 왼쪽사진은 그대로 두고 오른쪽 사진을 왼쪽 사진과 잘 맞게 변형을 할 것이다.
    그래서 코드의 전체 흐름은 오른쪽 사진을 image A로 두고 image A를 변환하는 흐름으로 진행된다.
  • [Line 8 ~ 22]
    앞선 글에서 다뤘던 내용이므로 설명은 생략하겠다.
  • [Line 25 ~ 27]
    cv2.findHomography( )를 통해 H를 구하는 과정이다.
    image A ( = 오른쪽 사진 )의 src ( = image A의 keypoint 좌표 )는
    image B ( = 왼쪽 사진 )의 dst ( = image B의 keypoint 좌표 )와 매칭되므로
    src가 dst로 변환될 수 있도록 해주는 변환 행렬 H를 cv2.findHomography( ) 함수가 구해준다.
    ※ cv2.RANSAC, status 에 대해서는 다른 글에서 정리하도록 하겠다.
  • [Line 30 ~ 35]
    변환 전 좌표들에 해당된다.
    의 모양에 맞춰주기 위한 과정이다.
  • [Line 37 ~ 40]
    H와 변환 전 좌표 before를 행렬곱 연산을 통해 변환 후 좌표 after를 얻는 과정이다.
    Line 38은 위에서 말한  을 나눠주는 과정에 해당하고,
    Line 39는 에서 와 만을 뽑는 과정이다.
  • [Line 43 ~ 53]
    최종적으로 변환 결과를 그려주는 과정이다.
    변환 후 좌표가 가능 범위를 벗어나는 경우에 에러가 날 수 있기 때문에 continue를 넣어줬다.
 

결과 이미지는 다음과 같다.

 

 

 왼쪽 이미지에 오른쪽 이미지가 맞춰서 변형된 걸 확인할 수 있다. 검정색 빗살이 쳐지는 이유는 아직 모르겠다. cv2.warpPerspective( ) 함수를 이용해 정합하면 저런 현상이 안 나타나는 걸 보면 코드에서 살짝 부족한 점이 있는 것 같다. 보간법 사용 여부에 따라 결과 차이가 생기는 건가?

 

 

 


! 우선 광고 시간 !

728x90
728x90

 


 

 

2019.10.08. 내용 추가

 

 검정색 빗살이 쳐지는 현상 해결에 대해 정리한 블로그를 우연히 읽게 됐다. 위와 같은 현상이 발생하는 이유는, before point가 after point로 변환되면서 after image의 모든 pixel 영역을 채워주지 못하기 때문이라고 한다. 그리고 이와 같은 문제는 보간(interpolation)을 통해 해결할 수 있다. max interpolation 방법이 가장 쉽게 구현할 수 있을 것 같아 구현해봤다.

h, w, _ = result.shape
kernel_size = 5
inpterpolation = np.zeros((h, w, 3), dtype = np.uint8)
for x in range(int((w / kernel_size))):
    for y in range((int(h / kernel_size))):
        roi = result[int(y * kernel_size) : int(y * kernel_size) + kernel_size,
                     int(x * kernel_size) : int(x * kernel_size) + kernel_size]

        roi_max_b = roi[:,:,0].max()
        roi_max_g = roi[:,:,1].max()
        roi_max_r = roi[:,:,2].max()
        
        for i in range(kernel_size):
            for j in range(kernel_size):
                if np.sum(roi[i, j] == [0,0,0]) == 3:
                    roi[i, j] = [roi_max_b, roi_max_g, roi_max_r]
        
        inpterpolation[int(y * kernel_size) : int(y * kernel_size) + kernel_size,
                      int(x * kernel_size) : int(x * kernel_size) + kernel_size] = roi
 
  • 각 채널 별로 max 값을 뽑아서 roi 이미지 중 [0, 0, 0]을 갖는 pixel에 대입하도록 구현했다.

 

 

 

 결과를 봐보면 검정색 빗살이 사라지긴 했는데 보간 결과가 그리 만족스럽진 않다. nearest 방식을 적용해보면 더 나으려나? 어떻게 구현하는 지 알게 되면 그때 또 추가로 정리해야겠다.

 

 

 


 

 

2020.07.10. 내용 추가

 

 마이클플리님이 댓글을 통해 남겨주신 조언을 토대로 알아본 결과 지금까지의 방식은 'Forward Mapping'이었고, 'Backward Mapping'을 통해 변환을 해주면 빗살 무늬 현상을 해결할 수 있다고 한다. 

 

[출처] https://www.researchgate.net/publication/323796200_Self-Supervised_Monocular_Image_Depth_Learning_and_Confidence_Estimation

 

 자세한 내용은 이 블로그를 참고하자. Backward Mapping은 다음과 같이 구현했다.

# backward mapping
before2 = []
for x in range(imageA.shape[1], imageA.shape[1] * 2):
    for y in range(imageA.shape[0]):
        point = [x, y, 1]
        before2.append(point)
before2 = np.array(before2).transpose()

Hinv = np.linalg.inv(H)

after2 = np.matmul(Hinv, before2)
after2 = after2 / after2[2, :]
after2 = after2[:2, :]
after2 = np.round(after2, 0).astype(np.int)

height, width, _ = imageA.shape
result2 = np.zeros((height, width * 2, 3), dtype = np.uint8)
for pt1, pt2 in zip(before2[:2, :].transpose(), after2.transpose()):
    if pt2[1] >= height or pt2[0] >= width:
        continue

    if np.sum(pt2 < 0) >= 1:
        continue
    
    result2[pt1[1], pt1[0]] = imageA[pt2[1], pt2[0]]
result2[0: height, 0 : width] = imageB
cv2.imshow('result2', result2)


# cv2.warpPerspective( )
result3 = cv2.warpPerspective(imageA, H, (imageA.shape[1] + imageB.shape[1], imageA.shape[0]))
result3[0 : imageA.shape[0], 0 : imageB.shape[1]] = imageB
cv2.imshow('result3', result3)
 
  • [Line 9]
    변환된 이미지의 좌표계를 기준으로 변환 전 이미지의 픽셀값을 가져오는 흐름이라
    Homography 행렬의 역행렬을 구해 변환에 사용하였다.
  • [Line 30 ~ 33]
    cv2.warpPespective( ) 함수를 통해 구한 결과 영상과 비교를 하기 위해 추가했다.

 

 

우선 Backward Mapping 결과 영상(result2)이다.

 

 

다음은 cv2.warpPespective( )를 사용했을 때의 결과 영상(result3)이다.

 

 

 Result2 영상이 보간을 사용했을 때 보다 훨씬 더 자연스러운 결과라는 것을 확인할 수 있다. 하지만 Result3에 비하면 살짝 부족한 느낌은 든다. Backward Mapping 단계에서 발생하는 소수형태의 좌표 결과를 어떻게 처리하느냐에 따라 영상의 퀄리티가 달라지는 것 같다. 나는 깊은 고민없이 단순히 반올림만 해준 것이라 아무래도 퀄리티가 상대적으로 낮은 것 같다. 나중에.. 좀 더 연구가 필요하다 생각이 들거나 우연히 또 알게된다면 시도해보독 하겠다.

 

 

 

[ 참고 사이트 ]

https://darkpgmr.tistory.com/79

https://b.mytears.org/2007/10/599/

https://wewinserv.tistory.com/89

https://yongchul-note.tistory.com/4

https://www.researchgate.net/publication/323796200_Self-Supervised_Monocular_Image_Depth_Learning_and_Confidence_Estimation

 

 

  • 마이클플리 2020.07.09 18:22

    좋은자료 감사드립니다.
    FindHomography 찾아보다 우연히 보고 답글 드리는데요.
    보간 하지 마시고 After Point좌표를 기준으로 역으로 Before에서 가져오시면 저런 현상 없어집니다.

    전에 회전 연산 했을때 나타났던 무늬랑 비슷한거 같은데요. 그때도 이미지 회전하면 회전전 영상의 픽셀을 회전후 로 맵핑 하면 int화 되면서 소수점이 날라가서 저렇게 되더라구요. 그래서 회전후 픽셀기준으로 회전전 픽셀을 가져오면 다 메워 집니다.

    • ballentain 2020.07.10 14:12 신고

      안녕하세요 마이클플리님!
      저도 정수화 과정에서 소수점이 날라가는 현상 때문에 발생하는 거라 생각은 하고 있었는데.. After Point 좌표 값을 기준으로 하는 방법까지는 생각하지 못 했었네요!! 좋은 조언 감사합니다ㅎㅎㅎ 시도해봤더니 보간보다는 확실히 괜찮은 결과가 나오네요!!

  • 애호박 2020.11.04 15:29

    안녕하세요. ballentain님

    올려주신 자료 정말 유익한 도움 되었습니다.
    파이썬에 대해 제대로 숙지하지 못한 채 이미지 정합을 하려다 보니 어려웠는데, 글을 보니 이해가 잘 되는 것 같습니다!

    궁금한 점이 있는데 두 이미지의 key point 점을 구하는 방식이 아닌 key point 점들의 좌표를 직접 입력해서 하는 방법도 있을까요??
    제가 정합하려는 이미지가 rgb 이미지와 열화상 이미지(흑백)이어서 그런지 key point 쌍이 아예 일치하지가 않더라구요ㅜㅜ

    • ballentain 2020.11.05 13:29 신고

      안녕하세요, 애호박 님!
      제 모자란 글이 도움이 되어서 기쁩니다ㅎㅎ

      RGB이미지와 열화상이미지 사이에서 매칭된 좌표들을 이미 알고 계신 상황이라면 cv2.getPerspectiveTransform( ) 함수를 사용하면 될 것 같습니다.

      ex)
      src = RGB 이미지에서 네 개 좌표
      dst = src에 대응하는 열화상 이미지에서 네 개 좌표
      M = cv2.getPerspectiveTransform(src, dst)
      warped_rgb_img = cv2.warpPerspective(RGB이미지, M, (width, height))

      # 변환된 RGB 이미지와 열화상 이미지 같이 표현
      result_img = cv2.addWeighted(warped_rgb_img , 0.5, 열화상이미지, 0.5, 0)

      cv2.getPerspectiveTransform( ) 함수는 제 블로그
      https://ballentain.tistory.com/40
      에 정리해뒀으니 한 번 읽어보시면 좋을 것 같습니다!

    • 애호박 2020.11.18 14:39

      정말 감사합니다! 다른 글들도 잘 읽어보겠습니다ㅎㅎ

  • 컴비 2020.12.10 14:06

    안녕하세요~ 너무 잘 봤어요! 혹시 첫번째 파트에서 h31x + h32y + 1 로 나누어야하는 이유를 알 수 있을까요..? 저렇게 하니까 진짜 되더라구요! 그런데 왜 되는지 이해는 안되네요ㅠ

    • ballentain 2020.12.11 11:18 신고

      안녕하세요, 컴비님! 저도 그 이유를 좀 알고싶네요.. 구글링 능력이 부족한 탓인지 아직 이유까지는 이해하기 못하고 그냥 외운 상황입니다ㅠㅡㅠ 나중에 알게되시면 댓글로 저도 좀 알려주세요!

    • 컴비 2020.12.12 02:56

      넵 알겠습니당ㅎㅎㅎ 다시 한번 감사드립니다!!

    • hsgo 2021.07.09 12:59

      homogeneous coordinates이기 때문입니다.

  • hsgo 2021.07.09 12:58

    bilinear interpolation과 같은 주변 픽셀 간의 보간을 사용하면 결과가 더 좋아질 것 같습니다.