해양기상을 위한 NetCDF 다루기

기본적인 python 개념을 알고 있다는 전제 하,
간단하게 사용하는 것들 위주로 기억하기 위해 최소한의 코드로 작성하는 글임을 참고 부탁드립니다.

---

NetCDF (Network Common Data Form) 데이터가 뭐에용?

해양과 대기 과학에서 널리 쓰인다.

• 다차원 데이터
• 시간과 공간에 따른 대기, 해양, 데이터들

배열 또는 행렬로 구성된 과학정보 형식

해양 수치 모델 결과는 공간정보(3차원)의 시간에 따른 값이 입력, 값 또한 여러 가지 변수를 포함. (x,y에 따른 h값. 면데이터)

또한 공간과 시간에 따른 파고, 파향, 파주기 등 부가 정보를 포함.

관측정보 일 때는 관측소명, 표준시, 수층, 장비명 등과 모델일 때는 시간, 공간 표준, 값 단위, fill_value 등등 포함되어야 할 정보들이 많음.

 

예로 시간과 공간 (x, y, t) + 그 값에 따른 z 까지면 4차원 데이터가 됨.

이런 시간에 따른 다수의 점 (공간, 명)의 값을 가지기에 NetCDF 형식을 활용하고 있음.

(대충 공간, 시간, 모든 정보 다차원 데이터라는 뜻)

NetCDF file structure

• Dimension (차원 : 위도, 경도, 시간의 크기, 갯수 등)
• Variables (실제 변수들의 값)
• attributes (파일에 대한 메타정보) 로 나누어져 있음.

확인하기 위해서는 → linux - ncdump, matlap - ncinfo

 

NetCDF Control

pip install netCDF4 seawater GSW xarray

conda install -c conda-forge basemap-data-hires

pip install netCDF4 seawater GSW xarray
conda install -c conda-forge basemap-data-hires

 

import numpy as np
import xarray as xr
import os
import matplotlib.pyplot as plt  # Corrected import for pyplot

# Set absolute path
work_path = '/DataProcessing_TRC/'
input_file = 'data/CTD_station_P1_NLEG01-1_-_Nansen_Legacy_Cruise_-_2021_Joint_Cruise_2-1.nc'
out_file = 'data/ctd_data_set.csv'

# Use path join to create full file paths
f_name = os.path.join(work_path, input_file)
f_out_name = os.path.join(work_path, out_file)

# Load the dataset
ds = xr.open_dataset(f_name)

# Extract temperature values
temperature = ds['TEMP'].values

# Convert selected variables to a DataFrame
df = ds[['PRES', 'TEMP', 'PSAL']].to_dataframe()

# Save the DataFrame to a CSV file
df.to_csv(f_out_name)

# Plotting
fig1 = plt.figure()

# Temperature plot
ax1 = fig1.add_subplot(121)
ax1.plot(ds['TEMP'], ds['PRES'], 'bo-', markersize=2, linewidth=0.5)
ax1.set_ylim(ax1.get_ylim()[::-1])  # Reverse the y-axis
ax1.set_xlabel('Temperature (C)')
ax1.set_ylabel('Depth (m)')
plt.grid(True)
ax1.xaxis.set_label_position('top')  # Move the label to the top
ax1.xaxis.set_ticks_position('top')   # Move the ticks to the top

# Salinity plot
ax2 = fig1.add_subplot(122)
ax2.plot(ds['PSAL'], ds['PRES'], 'ro-', markersize=2, linewidth=0.5)
ax2.xaxis.set_label_position('top')
ax2.xaxis.set_ticks_position('top')
ax2.set_ylim(ax2.get_ylim()[::-1])  # Reverse the y-axis
plt.grid(True)
ax2.set_xlabel('Salinity (PSU)')

# Show the plots
plt.show()

 

file 저장도 잘 된 것을 확인할 수 있다.

NetCDF 동해 평균해류도 다루기

format 확인하기 

# Set absolute path
work_path = '/mnt/c/Users/MML/sh/DataProcessing_TRC/'
input_file = 'data/KHOA_SCU_L4_Z004_Y01_U2023_EastSea.nc'

f_name = work_path + input_file
f_data = xr.open_dataset(f_name)

print(f_data)

u_comp = f_data.data_vars['u']
print(u_comp)

 

어떤 포맷인가 전반적으로 보고...

data variable들이 뭐가 있나 보기도 하고..

Realizing The East Sea NetCDF

''' File Processing '''
f_name = work_path + input_file
ds = xr.open_dataset(f_name)                    # Open the NetCDF dataset

# f_data = xr.open_dataset(f_name)
# print(f_data)
# u_comp = f_data.data_vars['u']
# print(u_comp)

lat = ds["lat"].values                          # 위도.values 로 numpy 변환
lon = ds["lon"].values                          # 경도.values 로 numpy 변환

# Transpose the u and v velocity components
u_velocity = np.transpose(ds["u"].values)       # East-west velocity
v_velocity = np.transpose(ds["v"].values)       # North-south velocity

• np.transpose()
  보통 해양 데이터는 시간, 깊이, 위도, 경도 등의 차원을 가지고 있어 속도 (u, v)의 배열이 시간과 깊이 차원을 먼저 포함한다면, 이를 전치하여 원하는 차원을 이용하기 위해 사용함!
  (대충 원하는 거 골라쓰려고 . 값을 가져오려고 뒤집었다는 뜻)

''' Convert to DataFrame for calculations '''
# Generate a grid for longitude and latitude
lon_grid, lat_grid = np.meshgrid(lon, lat)
# Create a DataFrame with longitude, latitude, and velocity components
data_frame = pd.DataFrame({
    "lon": lon_grid.ravel(),
    "lat": lat_grid.ravel(),
    "u_velocity": u_velocity.ravel(),           # Flatten to 1D
    "v_velocity": v_velocity.ravel()            # Flatten to 1D
})

• .ravel()
  원본이 아닌 view를 반환.
  다차원 -> 1차원으로 변환하여 Pandas DataFrame을 생성하고 있음.
  행 row 로 나열되어 쉽게 저장하고 이용.

''' Map plotting '''
fig, ax = plt.subplots(figsize=(10, 8))
m = Basemap(llcrnrlon=126,                      # Lower-left corner longitude
            llcrnrlat=33,                       # Lower-left corner latitude
            urcrnrlon=142,                      # Upper-right corner longitude
            urcrnrlat=48,                       # Upper-right corner latitude
            lon_0=134,                          # Central longitude
            lat_0=40,                           # Central latitude
            projection="merc",                  # Mercator projection
            resolution="h",                     # High resolution
            ax=ax)                              # Specify the axes object

# Fill continents and draw coastlines
m.fillcontinents(color="#D3D3D3", lake_color="lightblue")
m.drawcoastlines()

• Figure 전체그림, Axes 플롯의 축 = plt.subplots(figsize=(가로,세로))
• Basemap ( 왼아래경도, 왼아래위도, 오른위경도, 오른위위도, 중심경도, 중심위도, 투영법, 해상도, 축 )
• .fillcontinents()
• .drawcoastlines()

''' Calculate velocity of sea '''
# np.hypot (a, b) Calculate between a and b. 
# Calculate the horizontal speed using u and v components.
speeds = np.hypot(data_frame.u_velocity, data_frame.v_velocity)
 
# Plot arrows to represent velocity
X, Y = m(data_frame.lon, data_frame.lat)
quiver = m.quiver(X, Y, data_frame.u_velocity, data_frame.v_velocity, speeds, cmap="jet")

• np.hypot (a, b) :
  a, b 두 점 사이 거리 계산
  해수의 동서 u 방향, 남북 v 방향 속도를 이용해 수평 속도를 구함.
• m.quiver ( 시작점x, 시작점y,
  동서방향속도(x축속도)u,
  남북방향속도(y축속도)v,
  speed (벡터의 크기나 강도),
  cmap (컬러맵 지정인자)

''' Add remaining elements '''
# Add color bar
cbar = plt.colorbar(quiver, 
                    ax=ax,
                    orientation='horizontal',
                    shrink=0.4, 
                    aspect=12, 
                    pad=0.04)

# Draw meridians and parallels
m.drawmeridians(np.arange(126, 142, 3), labels=[0, 0, 0, 1])  # Meridians with label on top
m.drawparallels(np.arange(33, 49, 3), labels=[1, 0, 0, 0])    # Parallels with label on left

# Set title and show the plot
plt.title("Ocean Current Vector Map using Quiver", fontsize=14)
plt.show()

• plt.colorbar ( 객체, 축, orientation 방향(horizontal / vertical),
  shrink(막대 크기), aspect(종횡비. 10-20), pad(막대와 축 간격))
• 경선 .drawmeridians(범위, 레이블[상단, 우측, 하단, 좌측])
  위선 .drawparallels(범위,레이블[상단, 우측, 하단, 좌측])

오 .. 신기...

 

 

---

썸네일 이미지 참고 및 NetCDF 참고자료 ! (한 번은 읽어보면 좋겠다)

https://ncar.github.io/CESM-Tutorial/notebooks/resources/netcdf.html