'hough transform'에 해당되는 글 3건

  1. 2016.05.11 허프 변환 - 원(2/2)
  2. 2016.05.03 허프 변환 - 원(1/2)
  3. 2016.04.28 허프 변환 - 직선

허프 변환을 이용한 원의 검출 중에 반지름을 고정시키지 않은 방법에 대해서 정리한다.


기본 이론은 앞에서 정리했고 보완된 작업 위주로 정리한다. 


1) 반지름 범위 반영.

2) 스무딩.

3) 로컬 극대값 계산.

4) 근접 원들 합치기.

5) 결과.

6) 문제점.


1) 반지름 범위 반영

찾고자 하는 원의 크기를 범위로 지정할 수 있는데, 범위가 넓어질 수록 메모리와 계산량이 증가한다. 또한 반지름의 크기를 어느 정도 세밀하게 찾을 것인가를 지정할 수 있는데 이 역시 세밀하게 할 수록 부하가 심해질 것이다. 일단 예제는 반지름이 30~60 픽셀 원을 찾도록 하고 간격은 5 픽셀 정도로 했다. accumulator 공간이 반지름 변수의 추가로 3차원이 되었다.

rmin_idx = 30;

rmax_idx = 60;

r_step = 5;

accumulator = zeros( height+2*rmax_idx, width+2*rmax_idx, rmax_idx );


for r_idx=rmin_idx : r_step : rmax_idx

    for y_idx=1 : height

        for x_idx=1 : width

            if( edge_image(y_idx,x_idx) > 0 )

                for degValue=5 : 5 : 360

                    a_idx = round( x_idx-r_idx*cos(pi*degValue/180) );

                    b_idx = round( y_idx-r_idx*sin(pi*degValue/180) );

                    accumulator( b_idx+r_idx, a_idx+r_idx, r_idx ) = accumulator( b_idx+r_idx,                                                                                            a_idx+r_idx, r_idx )+1;

                end

            end

        end

    end

end

z축의 관점에서 보면 6군데에서 값이 크게 보인다. 

보이는 각도를 변경시키면 더 명확하게 확인할 수 있다. 반지름이 30, 45, 60정도에 주로 분포한다. 

각 반지름에서의 누적된 분포를 보면 다음과 같다. 







2) 스무딩

로컬 최대값을 찾기 위해서 스무딩처리를 한 후 임계값 이하는 0으로 처리한다. 

thresholdValue = 5.0;

for fig=0 : (rmax_idx-rmin_idx)/r_step

    figure(fig+1);

    

    smoothAccum(:,:,rmin_idx+fig*r_step) = imfilter( accumulator(:,:,rmin_idx+fig*r_step),                                                                       fspecial('disk',5), 'replicate' );

    

    for row=1:yd

        for col=1:xd

            if( smoothAccum(row,col,rmin_idx+fig*r_step) < thresholdValue )

      smoothAccum(row,col,rmin_idx+fig*r_step) = 0.0;

            end

        end

    end


    mesh( smoothAccum(:,:,rmin_idx+fig*r_step) );

    view(-15,70);

    title( strcat( 'radius : ', num2str(rmin_idx+fig*r_step) ), 'color', 'black' );

end


3) 로컬 극대값 계산

로컬 극대값을 찾는 함수를 사용해 극대값 위치를 찾는데 편평한( 모든 값이 일정한 ) 평면의  경우 모든 점을 로컬 극대값으로 계산하기 때문에 편평하지 않은 이미지에 대해서만 계산하도록 한다.

posCnt = 0;

posIdx = zeros(10,3); % ( x, y, radius )

for fig=0 : (rmax_idx-rmin_idx)/r_step

    posImg = imregionalmax( smoothAccum(:,:,rmin_idx+fig*r_step) );

    if( max(max(posImg)) == min(min(posImg)) )

        disp( 'no data' );        

    else

        for y=1 : yd

            for x=1 : xd

                if( posImg(y,x) > 0 )

                    posCnt = posCnt+1;

                    posIdx(posCnt,1) = x-(rmin_idx+fig*r_step);

                    posIdx(posCnt,2) = y-(rmin_idx+fig*r_step);

                    posIdx(posCnt,3) = rmin_idx+fig*r_step;

                end

            end

        end

    end    

end


4) 근접 원들 합치기.

같은 원을 다른 원으로 찾은 것 들을 하나로 합치는데 중심간의 거리가 최소 반지름보다 작은 것들은 하나라고 간주한다.( 중심이나 반지름이 약간 다른 원을 합치면서 정확한 원을 약간 벗어나게 찾을 수 있다. 이 부분은 개선의 여지가 있다. )

[ posCnt col ] = size( posIdx );  

cirIdx = posIdx;

cirIdx(:,4) = zeros();

cirCnt = 0;


for outer=1 : posCnt

if( cirIdx(outer,4) == 0 )

cirCnt = cirCnt+1;

cirIdx(outer,4) = cirCnt;

for inner=outer+1:posCnt

eucDist = sqrt( (cirIdx(outer,1)-cirIdx(inner,1))^2+(cirIdx(outer,2)-cirIdx(inner,2))^2 );

if( eucDist < rmin_idx )

cirIdx(inner,4) = cirIdx(outer,4);

end       

end

end

end


mergedMax = max( cirIdx(:,4) );

mergeIdx = 1;


mergeCnt = zeros( mergedMax, 1 );

xtotal = zeros( mergedMax, 1 );

ytotal = zeros( mergedMax, 1 );

rtotal = zeros( mergedMax, 1 );


for idx=1 : posCnt

xtotal( cirIdx(idx,4) ) = xtotal( cirIdx(idx,4) ) + cirIdx(idx,1);

ytotal( cirIdx(idx,4) ) = ytotal( cirIdx(idx,4) ) + cirIdx(idx,2);

rtotal( cirIdx(idx,4) ) = rtotal( cirIdx(idx,4) ) + cirIdx(idx,3);

mergeCnt( cirIdx(idx,4) ) = mergeCnt( cirIdx(idx,4) )+1;

end

for circle=1 : mergedMax

circleIndex(circle,1) = xtotal(circle)/mergeCnt(circle);

circleIndex(circle,2) = ytotal(circle)/mergeCnt(circle);

circleIndex(circle,3) = rtotal(circle)/mergeCnt(circle);

end



5) 결과



6) 문제점

① 이진화의 정밀도에 따라서 이후 작업에 영향 

    - 노이즈의 비율에 따라 끼치는 영향의 정량적인 분석 필요.

② 원주의 가려지는 비율과 노이즈 또는 threshold 값 결정에 주는 영향 정도 파악.

③ accumulator 공간에서 local maxima를 찾기위한 threshold 값을 결정하는 방법 고민.


'프로그래밍 > 영상처리' 카테고리의 다른 글

ubuntu octave 설치  (0) 2016.07.22
허프 변환 - 원(1/2)  (0) 2016.05.03
허프 변환 - 직선  (0) 2016.04.28
모달리티(Modality)  (2) 2016.04.13
Posted by 게 르 니 카

허프 변환을 이용한 원의 검출에 대해서 정리한다.


작업 순서는 다음과 같은데 설명은 기본 이론 요약 후 결과 이미지를 보며 설명한다.


1) 이론 이해.

2) 이미지 로딩.

3) 이진화.

4) 엣지 검출.

5) HT 공간 변환.

6) accumulator local maxima 찾기.

7) ( x, y )좌표에 원 표시.

8) 문제점.


1) 이론 이해.

일단 원을 수식으로 표현하면 다음과 같다.

원의 중심은 ( a, b ), 반지름은 r 인 원.


( x, y ) 평면을 ( a, b, r ) 공간으로 변환시키기 위해서 원의 식을 다른 방법으로 표현한다.

θ는 0 ≤ θ ≤ 2π .


( x, y ) 평면 상의 후보 픽셀을 ( a, b, r ) 공간으로 변환하는 작업을 아래와 같이 수행한다.


for r from r_min to r_max

for x from x_min to x_max

           for y form y_min to y_max

accumulator( a, b, r ) += 1

end

end

end


( a, b, r ) 공간의 값들 중에 지역 최대 값을 찾아 다시 ( x, y )평면으로 변환하면 원을 구할 수 있다.  


2) 이미지 로딩

일단 중간 크기 동전을 찾도록 반지름을 45으로 고정해서 계산한다. 이 후에 반지름 고정 없이

모든 원을 찾을 수 있도록 업그레이드 시킨다.

filename = 'coins3.jpg';

image = imread( filename );

3) 이진화 검출

컬러 이미지를 grayscale로 변환하고 Otsu's 방법을 이용해서 threshold 값을 찾고 이진화 후 'open' 모폴로지 연산을 적용하여 작은 조각들을 제거한다.

gray_image = rgb2gray( image );

level = graythresh( uint8(gray_image) );

binary_image = im2bw( uint8(gray_image), level );

binary_image = bwmorph( binary_image, 'open' );

4) 엣지 검출

sobel 필터를 사용해서 엣지 검출.

edge_image = edge( binary_image, 'sobel' );

5) HT 공간 변환.

윗 수식에 맞게 공간 변환 작업 수행.

수식 그대로 적용하면 싸인 함수 때문에 마이너스 인덱스를 갖을 수 있어서 에러가 발생한다.

최대 반지름 만큼 수평 이동시켜 인덱스를 구하고 나중에 값에서 최대 반지름 값을 빼준다.

accumulator = zeros( height+2*r_idx, width+2*r_idx );


for y_idx=1 : height

for x_idx=1 : width

if( edge_image(y_idx,x_idx) > 0 )

for degValue=5 : 5 : 360

a_idx = round( x_idx-r_idx*cos(pi*degValue/180) );

b_idx = round( y_idx-r_idx*sin(pi*degValue/180) );

accumulator( b_idx+r_idx, a_idx+r_idx ) = accumulator(                                                                                                b_idx+r_idx, a_idx+r_idx )+1;

end

end

end

end


6) accumulator local maxima 찾기.

아래 이미지의 가운데 상 하 원의 중심점 찾기.

thresholdValue = 5;

smoothAccum = imfilter( accumulator, fspecial('disk',5), 'replicate' );

[ height2 width2 depth2 ] = size( smoothAccum );


for y=1 : height2

for x=1 : width2

if( smoothAccum(y,x) < thresholdValue )

smoothAccum(y,x) = 0;

end

end

end

posImg = imregionalmax( smoothAccum );

posIdx = zeros(5,2);

posCnt = 0;


for y=1 : height2

for x=1 : width2

if( posImg(y,x) > 0 )

posCnt = posCnt+1;

posIdx(posCnt,1) = x-r_idx;

posIdx(posCnt,2) = y-r_idx;

end

end

end


7) (x,y)좌표에 원 표시.

원 표시.

NOP = 120;

THETA = linspace(0,2*pi,NOP);

RHO = ones(1,NOP)*r_idx;

[ origX, origY ] = pol2cart( THETA, RHO );


for circ=1 : posCnt

X = origX+posIdx(circ,1);

Y = origY+posIdx(circ,2);

plot(X,Y,'g', 'LineWidth', 2);

end


8) 문제점

① 이진화의 정밀도에 따라서 이후 작업에 영향 

    - 노이즈의 비율에 따라 끼치는 영향의 정량적인 분석 필요.

② 원주의 가려지는 비율과 노이즈 또는 threshold 값 결정에 주는 영향 정도 파악.

③ accumulator 공간에서 local maxima를 찾기위한 threshold 값을 결정하는 방법 고민.



'프로그래밍 > 영상처리' 카테고리의 다른 글

ubuntu octave 설치  (0) 2016.07.22
허프 변환 - 원(2/2)  (0) 2016.05.11
허프 변환 - 직선  (0) 2016.04.28
모달리티(Modality)  (2) 2016.04.13
Posted by 게 르 니 카

허프 변환은 이미지 속에서 주요 특징 요소 들을 - 직선, 원, 타원, 등등 - 찾는 방법이다.

1959년 Paul Hough가 'bubble chamber' 사진의 기계적 분석을 위해 고안했으며,

1962년 "Method and Means for Recognizing Complex Patterns" 이란 타이틀로 미국 특허를 

획득했다. 더 자세한 이력이나 설명은 위키에 자세히 나와 있다.


우선 기본 원리부터 이해하자.


직선을 수학적으로 표현하는 방법에는 여러 가지가 있다.


( x, y ) 평면에서 기울기 m, 절편 n인 직선은

y = m*x + n.


( x1, y1 ), ( x2, y2 ) 두 점을 지나는 직선은

y - y1 = ( ( y2-y1 ) / ( x2-x1 ) ) *( x - x1 ).


x 절편이 a, y 절편이 b인 직선은 아래와 같다.

아래 그림을 보자.

 


빨간 직선을 절편이 아닌 각도와 거리로 표현하려면 어떻게 해야할까?

빨간 직선의 x 축과 만나는 점을 x1, y 축과 만나는 점을 y1이라 하면  

( x, y ) 평면 상의 한 점 ( x2, y2 )를 위 직선의 방정식에 대입하면

r = x2cosθ +y2sinθ 이 되어 θ 값에 따라 r의 값이 변하는 곡선이 된다.


x 축은 θ가 0이고 -x 축은 180도가 되어 θ가 0 <= θ < pi/2 이 구간이면 ( x2, y2 ) 점을 지나는 

( x, y ) 평면상의 모든 직선을 표현할 수 있다. 


아래 그림을 보면

첫 화면은 (2,2)점을 지나는 기울기가 다른 직선들을 표시한 것이다. 초록색 직선의 경우 θ는 30도 정도 될 것이다.( 원점을 지나며 초록색 직선과 수직인 직선의 x축과의 각도 ), 빨간색 직선은 60도.


두번째 화면은 각 직선들의 ( θ, r ) 평면상의 위치를 표현한 것으로 r = 2cosθ +2sinθ 인 직선에서 각도가 30도로 고정되면 r의 길이도 결정되어 한 점이 된다.


세번째는 각도를 0도 부터 180도까지 1도씩 증가시킨 것으로 sine 곡선이 된다.


( x, y ) 평면 상의 5개 점을(왼쪽 그림) ( θ, r ) 평면에 표현하면(오른쪽 그림) 0도에서 180도로 변함에 따라 r의 값이 변하여 삼각함수 곡선이 5개가 그려진다. 5개의 점이 일직선 상에 있다면 ( θ, r ) 평면의 5개 곡선은 한 점에서 만나게 된다. 위 예의 경우 만나는 지점의 위치는 ( 45, 2.828 ) 이다.  r = xcosθ +ysinθ  직선에 θ=45, r=2.828을 대입하면 2.828 = xcos(45)+ysin(45) 는  y = -x + 4 의 직선을 얻을 수 있다.

'프로그래밍 > 영상처리' 카테고리의 다른 글

ubuntu octave 설치  (0) 2016.07.22
허프 변환 - 원(2/2)  (0) 2016.05.11
허프 변환 - 원(1/2)  (0) 2016.05.03
모달리티(Modality)  (2) 2016.04.13
Posted by 게 르 니 카
이전버튼 1 이전버튼

블로그 이미지
게 르 니 카

공지사항

Yesterday
Today
Total

달력

 « |  » 2024.5
1 2 3 4
5 6 7 8 9 10 11
12 13 14 15 16 17 18
19 20 21 22 23 24 25
26 27 28 29 30 31

최근에 올라온 글

최근에 달린 댓글

최근에 받은 트랙백

글 보관함