작업5 : Fortran으로 합성 풍향 풍속 구하기1

 

이번에는 합성 풍속과 풍향을 구해보겠습니다

 

평균 풍속 풍향이라고만 하면 하루동안 측정한 자료를 모아서 24시간으로 나누면 되지 않냐 생각하실 수도 있습니다

 

풍속은 1시간마다 측정하여 24로 나누는 산술평균 방법을 사용해도 괜찮지만

 

방향을 포함하고 있는 풍향은 얘기가 좀 다릅니다

 

벡터 형태이기 때문에 방향이 정말 중요합니다

---

바람의 성분)

바람 벡터를 동서 성분과 남북 성분으로 나누어서 주로 생각합니다

 

쉽게 표현하기 위해 동서 성분은 u성분, 남북 성분을 v성분이라고 부릅니다

 

합성 풍향 풍속을 구하기 위해서는 평균 풍향 풍속을 먼저 구해야합니다

 

원래는 이런 공식을 사용해서 구한 다음 해야하지만 조금 더 쉬운 공식을 사용하겠습니다

 

그전에 기상학적 바람 벡터와 수학적 바람벡터부터 알아보겠습니다

---

기상학적 바람 벡터와 수학적 바람 벡터)

.

자 여기 기상학적으로 표현된 바람 벡터 그림이 있습니다

 

기상학적으로는 풍향을 나타낼 때 불어오는 곳을 기준으로 표현합니다

 

s 방향에서 불어가면 남풍, w 방향에서 불어가면 서풍 이렇게 말이죠

 

북쪽을 0도로 잡고 시계방향으로 각도가 증가합니다

 

북쪽이 0도, 동쪽이 90도, 남쪽이 180도, 서쪽이 270도

 

하지만 이걸 수학적인 벡터로 표현하면 많이 달라집니다

 

수학적으로 바람을 표현하면 불어가는 방향을 나타내는 것입니다

 

사방위가 반대인 이유는 위에서 아래로 부는 바람이 북풍이기 때문에

 

북풍 벡터다 라는 뜻을 보여주기 위한 겁니다

 

각도 표시도 조금 다릅니다

 

기상학적 바람에서는 북풍이 0도 였지만 수학적으로는 사분면 개념을 사용하기 때문에

 

서쪽을 기준으로 반시계 방향으로

 

서쪽이 0도, 남쪽이 90도, 동쪽이 180도, 북쪽이 270도 이렇게 표현을 합니다

 

지금까지 이걸 왜 설명했냐 하면 포트란에 바람의 각도를 입력시킬때는 수학적인 바람표현이 필요하지만

 

기상청이 제공하는 자료는 기상학적인 바람 표현이기 때문입니다

 

동쪽을 예로 들어서 설명해보면 기상학적인 바람표기는 0도 이지만 수학적 표기는 270도인 것이죠

 

이걸로 봤을 때 각도 변환을 270 - 기상학적 바람을 해주어야 한 다는 것을 알 수 있습니다

 

이렇게 변환된 각도를 가지고

 

위의 공식을 이용해서 u성분과 v성분을 구할 수 있는것이죠

 

이번 작업은 개념 설명도 길고 코드도 쪼금 길어서 포스팅을 2개로 나눠서 하겠습니다

 

너무 길면 읽기 힘드니까요

 

첫번째에서 코드에 각도가 왜 저렇게 들어갔는지 이해하시고

 

두번째로 넘어오시면 좋을 것 같습니다

 

혹시 코드를 개인적으로 해석해보실 분들을 위해 코드는 마지막에 적어두겠습니다

 

        REAL, DIMENSION(24) :: ws, wd

        REAL, DIMENSION(24) :: u, v, cpwd
        REAL, parameter :: pi = 4. * atan(1.)

        REAL :: ap1, ap2, ap3, ap4

        CHARACTER(2), parameter :: ug=char(int(Z'C2'))//char(int(Z'B0'))

        open(1, file = 'time.dat')
        open(2, file = 'time2.dat')

        write(2, 60) 'u component', 'v compnent', 'compose component'
60      format(a, 3x, a, 3x, a)

        do i = 1, 24
        read(1, 10) ws(i), wd(i)
10      format(6x, f3.1, 1x, f3.0)


        u(i) = ws(i) * cos((270. - wd(i)) * (pi/180))
        v(i) = ws(i) * sin((270. - wd(i)) * (pi/180))

        cpwd(i) = atan2(v(i), u(i))

        ap1 = ap1 + u(i)
        ap2 = ap2 + v(i)
        ap3 = ap3 + cpwd(i)

        write(2, 20) u(i), v(i), cpwd(i)

20      format(f10.4, 2x, f10.4, 5x, f10.4)

        end do

        rp1 = ap1/24.
        rp2 = ap2/24.
        rp3 = (270. - ((atan2(rp2, rp1)) * (180/pi))) - 360.

        write(2, 30) "compose u direction componet :", rp1
30      format(a,f10.4, 1x, 'm/s')

        write(2, 40) "compose v direction component :", rp2
40      format(a, f10.4, 1x, 'm/s')

        write(2, 50) "compose wind speed :", sqrt(rp1**2 + rp2**2)
50      format(a,f10.4, 1x, 'm/s')

        write(2, 70) 'compose wind direction', rp3, ug
70      format(a, f8.2, 1x, a)
      
        end