drag1.nb
0.02MB

공기 저항이 없을 때 물체를 $v_0$ 속력으로 위로 던지면 최고점에 올라가는데 걸리는 시간과 다시 내려오는데 걸리는 시간은 동일하게 $t_{ff} = v_0/g$로 주어진다. 공기 저항이 있는 경우는 어떻게 될지 구체적으로 계산해보자.

반지름이 $R$인 공 모양의 물체가 속력의 제곱에 비례하는 공기 저항(끌림힘) $D= \frac{1}{2} C\rho_{air} A v^2 \approx  0.2 \rho_{air} \pi R^2 v^2$을 받을 때, 올라가는 동안 운동 방정식은 (물체가 받는 공기의 부력도 고려해야 하지만 여기서는 무시한다. 부력은 $g$을 약간 줄이는 효과를 만든다)

$$  m \frac{dv}{dt} = -mg -\ 0.2   \rho_{air} \pi R^2  v^2, $$

$$ \rightarrow  \quad \frac{dv}{dt} = -g (1 + \gamma^2 v^2) , \quad\quad  \gamma^2 =  0.2 \frac{\pi R^2 \rho_{air} }{mg} .$$

$\gamma$는 내려오는 과정에서 terminal speed의 역수를 의미한다. 시간에 대해 적분을 해서 속도를 구하면,

$$ v(t) = \frac{1}{\gamma} \tan \Big( -\gamma gt + \arctan(\gamma v_0) \Big)$$

최고점에 도달하는데 걸리는 시간은 $v(t_{up})=0$ 에서

$$ t_{up} = \frac{1}{\gamma g} \arctan(\gamma v_0) = \frac{v_0}{g} \Big(1 - \frac{(\gamma v_0)^2}{3}+....\Big) $$

로 주어지므로 공기저항이 없을 때보다 더 짧다. 최고점의 높이는

$$ h_{max} = \int_0^{t_{up}} v(t) dt = \frac{1}{\gamma^2 g} \ln \sqrt{ 1 + (\gamma v_0)^2 } =\frac{v_0^2}{2g} \Big( 1- \frac{(\gamma v_0)^2}{2 }+...\Big)$$

로 주어지므로 역시 공기 저항이 없을 때보다 낮다.

 

다시 내려오는 과정에서 중력과 끌림힘이 반대방향이므로 운동 방정식은

$$ \frac{dv}{dt}= -g (1 -\gamma^2 v^2)$$

로 주어지고, 이를 적분하면 (출발 시간을 $t=0$으로)

$$ v(t) = -\frac{1}{\gamma} \tanh (\gamma g t)$$

이를 다시 적분하면 낙하시간에 따른 높이를 얻을 수 있다: $h(0)=h_{max}$

$$ h(t) = \frac{1}{\gamma^2 g} \ln \frac{ \sqrt{ 1+ (\gamma v_0)^2 } }{ \cosh( \gamma gt)}$$

따라서 바닥에 떨어지는데 걸리는 시간 $t_{dn}$은: $ h(t_{dn}) = 0$

$$t_{dn} = \frac{1}{\gamma g}\text{arccosh}\sqrt{1+(\gamma v_0)^2} =\frac{1}{\gamma g} \ln \left(\gamma v_0 + \sqrt{1+ (\gamma v_0)^2 } \right) = \frac{v_0}{g}\Big( 1  -  \frac{ (\gamma v_0)^2 }{6}+...\Big)$$

 

두 시간을 비교해보면 위로 던져진 물체가 최고점에 올라가는데 걸리는 시간보다 다시 내려오는데 걸리는 시간이 더 길다는 것을 볼 수 있다:

$$ {t_{dn} - t_{up}} \approx  \frac{(\gamma v_0)^2}{6} t_{ff}$$

 

$v_0 =10\text{m/s}$로 ($\rightarrow t_{ff} = 1.02\text{s}$, $v_{terminal} \approx 36.6\text{m/s}$) 야구공을 공중으로 던지는 경우를 예로 들면, 차이는 대략 0.013초 정도로 계산된다. 던지는 속력이 종단속력에 가깝거나 더 크면 근사식을 사용할 수 없고 정확한 계산식을 이용해야 한다.

 

$v_0 = 40\text{m/s}$, $1/\gamma=36.6 \text{m/s}$ 일 때,

728x90
Posted by helloktk
,

이미지 처리 과정에서 미분에 해당하는 그래디언트 필드(gradient field: $g_x$, $g_y$ )를 이용하면 이미지 상의 특징인 corner, edge, ridge 등의 정보를 쉽게 얻을 수 있다. 이미지의 한 지점이 이러한 특징을 가지는 특징점이 되기 위해서는 그래디언트 필드의 값이 그 점 주변에서 (3x3나 5x5정도 크기의 window) 일정한 패턴을 유지해야 한다. 이 패턴을 찾기 위해서 그래디언트 필드에 PCA를 적용해보자. 수평과 수직방향의 그래디언트 field인 $g_x$와 $g_y$ 사이의 covariance 행렬은 다음 식으로 정의된다:

$$ \Sigma = \left [ \begin {array}{cc} < g_x^2 > & < g_x g_y > \\    <g_x g_y> & <g_y^2 > \end {array}\right] =\left [ \begin {array}{cc} s_{xx} & s_{xy} \\ s_{xy} & s_{yy}\end {array}\right];$$

$<...> = \int_{W}(...) dxdy$는 픽셀 윈도에 대한 적분을 의미한다. $\Sigma$의 eigenvalue는 항상 음이 아닌 값을 갖게 되는데 (matrix 자체가 positive semi-definitive), 두 eigenvalue이 $λ_1$, $λ_2$면

$$λ_1 + λ_2 = s_{xx} + s_{yy} \ge 0, \quad \quad    λ_1  λ_2 = s_{xx} s_{yy} - s_{xy}^2 \ge0 $$

을 만족한다 (완전히 상수 이미지를 배제하면 0인 경우는 없다). eigenvalue $λ_1$, $λ_2$는 principal axis 방향으로 그래디언트 필드의 변동(분산)의 크기를 의미한다. edge나 ridge의 경우는 그 점 주변에서 잘 정의된 방향성을 가져야 하고, corner의 경우는 방향성이 없어야 한다. edge나 ridge처럼 일방향성의 그래디언트 특성을 갖거나 corner처럼 방향성이 없는 특성을 서로 구별할 수 있는 measure가 필요한데, $λ_1$과 $λ_2$를 이용하면 차원이 없는 measure을 만들 수 있다. 가장 간단한 차원이 없는 측도(dimensionless measure)는  eigenvalue의 기하평균과 산술평균의 비를 비교하는 것이다.

$$ Q = \frac { {λ_{1} λ_{2}} }{ \left( \frac {λ_{1}+λ_{2}}{2} \right)^2} = 4\frac { s_{xx} s_{yy} - s_{xy}^2}{(s_{xx} + s_{yy})^2};$$

기하평균은 산술평균보다도 항상 작으므로

$$ 0 \le Q \le 1 $$

의 범위를 갖는다. 그리고 $Q$의 complement로

$$P = 1-Q = \frac{(s_{xx}-s_{yy})^2 + 4 s_{xy}^2}{(s_{xx}+s_{yy})^2};$$를 정의할 수 있는 데 $0 \le P \le 1$이다. $Q$와 $P$의 의미는 무엇인가? 자세히 증명을 할 수 있지만 간단히 살펴보면 한 지점에서 $Q \rightarrow 1$이려면 $λ_{1} \approx λ_{2}$이어야 하고, 이는 두 주축이 동등하다는 의미이므로 그 점에서는 방향성이 없는 코너의 특성을 갖게 된다. 반대로 $Q \rightarrow 0$이면 강한 방향성을 갖게 되어서 edge(ridge) 특성을 갖게 된다.

 

실제적인 응용으로는 지문 인식에서 지문 영역을 알아내거나 (이 경우는 상당이 큰 윈도를 사용해야 한다) 또는 이미지 텍스쳐 특성을 파악하기 위해서는 이미지를 작은 블록으로 나누고 그 블록 내의 미분 연산자의 균일성을 파악할 필요가 있는데 이 차원이 없는 측도는 이미지의 상태에 상관없이 좋은 기준을 주게 된다.

 

참고 논문:

Image field categorization and edge/corner detection from gradient covariance
Ando, S.

Pattern Analysis and Machine Intelligence, IEEE Transactions on
Volume 22, Issue 2, Feb 2000 Page(s):179 - 190

 

** 네이버 블로그 이전;

728x90

'Image Recognition > Fundamental' 카테고리의 다른 글

Is Power of 2  (0) 2021.02.12
Flood-Fill and Connected Component Labeling  (2) 2021.02.10
점증적인 cosine/sine 값 계산  (0) 2020.12.28
Fast Float Sqrt  (0) 2020.12.27
2차원 Gaussian 분포 생성  (0) 2020.12.10
Posted by helloktk
,

기하 알고리즘을 특히 폴리곤 관련 알고리즘을 테스트하기 위해서는 폴리곤 데이터를 만들어야 한다. 여기서 2차원의 점들을 이용하여 간단히 simple polygon을 만드는 방법을 소개한다.

 

단순 폴리곤을 만들기 위해서는 점들을 정렬을 해야 한다. 각각의 점들을 x 축에 프로젝션을 시키면 x 축에 대해서 정렬이 된다(같은 x 값인 경우에는 y의 크기대로 정렬). 따라서 x 값으로 정렬된 점들을 쭉 이어나가면 하나의 poly-line을 만들 수 있다. poly-line의 끝을 잇는다고 해서 단순 폴리곤이 항상 만들어지지는 않는다. 이를 해결하기 위해서는 점들을  한 직선을 기준으로 위-아래로 분리하고, 직선의 윗부분은 x 값이 큰 순서대로 정렬을 하고, 아래 부분은 x 값이 작은 순서대로 정열을 하면 폴리 라인이 아닌 단순 폴리곤을 얻을 수 있다. 직선은 어떻게 잡을까? 이것은 가장 작은 x 값을 갖는 점과 가장 큰 x 값을 같은 점을 잇는 직선을 생각하면 편리하다.

// cross(BC, BA)
// A->B-C: 반시계(C가 AB 왼편: > 0), 시계(C가 AB오른편: < 0), 일직선(0);
int CCW(CPoint A, CPoint B, CPoint C) {
    C.x -= B.x; C.y -= B.y; 
    A.x -= B.x; A.y -= B.y;
    return C.x * A.y - C.y * A.x;
}
int comp(const void *A, const void *B) { //x->y의 크기 순서대로 정렬;
    int ax = ((CPoint *)A)->x - ((POINT *)B)->x ;
    if (ax > 0) return  1;
    if (ax < 0) return -1;
    int ay = ((CPoint *)A)->y - ((CPoint *)B)->y ;
    if (ay > 0) return  1;
    if (ay < 0) return -1;
    return 0;
}
void MakeSimplePolygon(std::vector<CPoint>& pts, std::vector<CPoint>& spoly) {
    if (pts.size() < 1) return;
    int maxid = 0, minid = 0;
    for (int i = pts.size(); i-->1;) {
        if (pts[i].x > pts[maxid].x) maxid = i;
        if (pts[i].x < pts[minid].x) minid = i;
    }
    std::vector<CPoint> LP, RP;
    for (int i = 0; i < pts.size(); i++) {
        if (i == maxid || i == minid) continue;
        // 기준선의 왼편(LP)/오른편(RP)에 놓인 점 분리;
        if (CCW(pts[minid], pts[maxid], pts[i]) >= 0) LP.push_back(pts[i]); 
        else RP.push_back(pts[i]);
    }
    if (LP.size() > 0) 
        qsort(&LP[0], LP.size(), sizeof(CPoint), comp);
    if (RP.size() > 0)    
        qsort(&RP[0], RP.size(), sizeof(CPoint), comp);
	
    spoly.clear();
    spoly.push_back(pts[minid]);
    for (int i = 0; i < LP.size(); i++) spoly.push_back(LP[i]);
    spoly.push_back(pts[maxid]);
    for (int i = RP.size(); i-- > 0;) spoly.push_back(RP[i]);
    // spoly = clockwise simple polygon;
};

다른 방법으로는 위와 같은 직선을 만들고, 아랫부분의 점들의 convex hull을 구성하고, 나머지 점들을 마찬가지의 방법으로 정렬을 하여 폴리곤을 만들거나, 아니면 한 점을 잡아서(맨 아래 오른쪽) 그 점을 기준으로 해서 나머지 점들을 각도에 대해서 정렬을 한 다음 그 정렬된 순서대로 폴리곤을 구성하면 된다.

기준점에 대한 각도를 정렬하여서 만든 예(동일한 각도가 생기는 경우에는 기준점에서의 거리로 비교)

**네이버 블로그에서 이전;

 
728x90

'Computational Geometry' 카테고리의 다른 글

Minimum Bounding Rectangle  (3) 2021.03.15
Minimum Enclosing Circle  (0) 2021.03.01
단순 다각형의 Convex hull  (0) 2021.01.24
단순 다각형의 무게중심  (0) 2021.01.24
단순 다각형의 면적(2D)  (0) 2021.01.23
Posted by helloktk
,

평면 상에 주어진 점들의 convex hull을 구하는 알고리즘 중의 하나인 Graham scan에서는 먼저 주어진 점들을 한 점을 기준으로 각도로 정렬하는 과정이 필요했다. 그러면 점들이 순차적으로 연결된 단순 다각형에서는 sorting과정이 없이 Graham scan 알고리즘을 바로 적용이 가능한 것일까? 이에 대한 counter example은 쉽게  찾을 수 있다. 단순 폴리곤에 대해서도 항상 각도에 대한 정렬 과정이 필요한 것인가? 답은 아니다. 정렬 과정이 없이도 단순 폴리곤에 대해서는 쉽게 convex hull을 구할 수 있는 알고리즘이 존재한다. 정렬 과정이 없이 단순 폴리곤의 convex hull을 찾는 알고리즘에 대해서 알아보자.

Melkman Algorithm;

우선 폴리곤의 이웃하는 세 꼭짓점을 잡아서, 반시계 방향의 삼각형을 구성한다. 이 삼각형을 deque에 넣는다 (bottom = top). 폴리곤을 순환하면서 새 점이 들어올 때 이미 만들어진 convex hull의 내부점이 아니면 이 점이 포함되도록 convex hull을 업데이트한다: Graham scan 알고리즘의 scanning 과정을 bottom을 기준으로 반시계 방향으로 convexity를 만족시킬 때까지 bottom을 제거하고, top을 기준으로는 시계방향으로 convexity를 만족시킬 때까지 top을 제거한다. 이 과정이 끝나면 새 점을 deque에 추가한다. 이 과정을 나머지 모든 점들에 대해서도 진행한다.

int LEFTSIDE(CPoint A, CPoint B, CPoint C){ // cross(AB, AC) > 0: C가 AB 대해서 왼쪽
	return ((B.x - A.x) * (C.y - A.y) - (B.y - A.y) * (C.x - A.x)) > 0;
}
int convexhull_spolygon(std::vector<CPoint>& pts, std::vector<CPoint>& hull) {
    int N = pts.size();
    if (N < 3) return 0;
    std::vector<CPoint> Q(2 * N + 1);  // deque;
    int bot = N - 2, top = bot + 3;
    // |0|1|......|N-2|N-1|N+0|N+1|N+2|..............|2N|;
    //              bot    top    
    Q[bot] = Q[top] = pts[2];
    if (LEFTSIDE(pts[0], pts[1], pts[2])) { //2->0->1->2; Q에 ccw로 입력;
        Q[bot + 1] = pts[0]; Q[bot + 2] = pts[1];
    } else {
        Q[bot + 1] = pts[1]; Q[bot + 2] = pts[0];
    }
    int i = 2; // last touched index;
    while (++i < N) {
        // 기존의 convex_hull에 들어있으면 제외.
        if (LEFTSIDE(Q[bot + 0], Q[bot + 1], pts[i]) && 
            LEFTSIDE(Q[top - 1], Q[top + 0], pts[i]))
            continue;
        // bottom에 대해서 ccw 방향으로 체크(bot 증가 방향)
        // pts[i]가 (bot)->(bot+1)라인의 오른쪽에 있으면, 기존의 bottom을 제외;
        while (!LEFTSIDE(Q[bot], Q[bot + 1], pts[i]) && (bot < top)) bot++;
        Q[--bot] = pts[i];
        // 추가점에서 top->top-1의 ccw 방향으로 convexity 체크하여 만족하지 않은 경우
        // top을 감소
        while (!LEFTSIDE(Q[top - 1], Q[top], pts[i]) && (top > bot)) top-- ;
        Q[++top] = pts[i];
    }
    hull.resize(top - bot);
    for (int i = hull.size() - 1; i >= 0; --i) hull[i] = Q[i + bot];
    return hull.size();
};

**네이버 블로그 이전;

 

 

728x90

'Computational Geometry' 카테고리의 다른 글

Minimum Enclosing Circle  (0) 2021.03.01
Creating Simple Polygons  (0) 2021.01.25
단순 다각형의 무게중심  (0) 2021.01.24
단순 다각형의 면적(2D)  (0) 2021.01.23
삼각형 외접원의 Inclusion Test  (0) 2020.12.30
Posted by helloktk
,

단순 다각형(simple polygon)의 무게중심(center of gravity or center of mass)은 다각형을 균일한 밀도의 판으로 생각했을 때 판의 무게중심과 같다. 가장 단순한 다각형인 삼각형의 무게중심은 세 꼭짓점의 산술평균으로 표현된다.

$$ \text{CoG} = \frac{1}{3} ({{\bf P} + {\bf Q} + {\bf R}}).$$

증명: 삼각형의 한 변 PQ에 나란한 띠로 삼각형을 분할하자. 그러면 각 띠의 무게중심은 띠의 기하학적 중심이므로 꼭지점 R와 변 PQ의 중심을 연결한 선분 RA상에 있어야 한다. 그리고 띠의 무게중심은 그 점에 모든 질량이 뭉친 것으로 생각할 수 있으므로 전체 삼각형의 무게중심은 선분 RA상의 어느 지점에 있어야 한다. 마찬가지 논리로  삼각형을 선분 QR에 나란한 띠로 분할나누면 이 띠들의 무게중심은 선분 PB상에 있음을 알 수 있다. 따라서 삼각형의 무게중심도 PB상의 한점이어야 하므로 선분 PB와 선분 RA의 교점이 무게중심이 된다. 

다각형은 삼각형으로 분할되므로 이 분할된 삼각형의 무게중심을 이용하면 쉽게 계산할 수 있다. 분할된  삼각형의 무게중심을 면적으로 가중치를 준 평균값이 다각형의 무게중심이 된다. 

 

실제 계산에서는 다각형을 삼각 분할하지 않고도 간단한 방법에 의해서 무게중심을 구할 수 있다. 원점과 다각형의 각 변의 꼭짓점을 이용해서 삼각형들을 구성하면 원래의 다각형을 겹치게(원점이 내부에 있으면 겹침이 없다) 분할할 수 있다. 분할된 삼각형으로 무게중심을 구할 때 겹치는 영역의 기여를 제거해야 한다. 그런데 다각형 밖의 영역을 분할하는 삼각형은 다각형 내부를 분할하는 삼각형과는 다른 orientation을 가지게 된다. 삼각형의 면적은 한 꼭짓점을 공유하는 두 변의 외적에 비례하므로, 반대의 orientation을 갖는 삼각형은 자동으로 반대 부호의 면적을 가지게 된다. 따라서 분할된 삼각형 면적 가중치를 외적으로 주어서 무게중심을 구하면 겹치는 영역이 자동으로 상쇄되는 효과를 얻을 수 있다.

노란색 영역에 있는 분할 삼각형은 음의 면적

$$\begin{align} \text{CoG} &= \frac{1}{\text{다각형 면적}} \sum (\text{삼각형 CoG})  (\text{면적})  \\ &= \frac{1}{\text{다각형 면적}} \sum \frac{1}{3} \left( {\bf V}_i + {\bf V}_{i+1} + {\bf O}\right ) \frac{ {\bf V}_{i} \times {\bf V}_{i+1} }{2} \\ &= \frac{1}{3}\frac{1}{   \text{다각형 면적} }\sum ( {\bf V}_{i} + {\bf V}_{i+1}) \frac{{\bf V}_{i} \times {\bf V}_{i+1}}{2} \end{align}$$ 

다각형의 면적($=\sum \frac{1}{2}({\bf V}_i \times {\bf V}_{i+1})$)을 구할 때 삼각형과 동일하게 orientation에 따라 부호를 포함하도록 설정하면 다각형의 면적 부호가 삼각형의 면적 부호로 상쇄되므로 다각형의 orientation에 무관하게 성립하는 공식이 된다.

CPoint polygon_centroid(CPoint V[], int N) {
    double cx = 0, cy = 0, area2 = 0;
    for(int i = 0, j = N - 1; i < N; j = i++) {
        double tri_area2 = V[i].x * V[j].y - V[i].y * V[j].x; // area * 2;
        cx += (V[i].x + V[j].x) * tri_area2;
        cy += (V[i].x + V[j].x) * tri_area2;
        area2 += tri_area2;                                   //total area * 2
    }
    cx /= 3 * area2;
    cy /= 3 * area2;
    return CPoint(int(cx + 0.5), int(cy + 0.5))
};

** 네이버 블로그에서 수정 이전;

728x90

'Computational Geometry' 카테고리의 다른 글

Creating Simple Polygons  (0) 2021.01.25
단순 다각형의 Convex hull  (0) 2021.01.24
단순 다각형의 면적(2D)  (0) 2021.01.23
삼각형 외접원의 Inclusion Test  (0) 2020.12.30
Point in Polygon  (2) 2020.12.14
Posted by helloktk
,