전세계 기온 시각화

지구 위에 데이터의 흐름을 올려본 일들은 있었는데 그건 C++에서 했던 작업이었다. GPU의 VRAM 용량이 충분하니 지구 전체를 한 판으로 만들어 돌렸었고, 미리 VRAM에 복사해놓은 데이터의 주소들을 여러시점에 따라 다르게 참조해서 그리는 방식이었다.

 

비슷해보이는 표현이라도 웹으로 구현하려고 하면 여러가지 제약이 따른다. 가장 첫 번째 제약은 데이터의 용량이다. 특정 시점에서화면에 그리는 데이터가 몇 MB를 넘어가버리게 되면, 아니, 1MB 정도라도 전송속도로 인해 끊김없이 보기에는 제약이 따른다. 물론 이런 상황에 대한 해법으로는 '타일맵'이라는 기술이 오래전에 나와 있다. 단, 이번에 한 '전 세계 기온 시각화' 작업에서는 시간 변화에 따라 3시간 간격의 원본 데이터를 자엽스럽게 이어서 보여주는 것을 목표로 했다.

 

 

데이터

 

전 세계의 기온 데이터로는 향후 일주일간의 예측을 포함한 GFS(Global Forecast System)의 예측 데이터가 있다.

ERDDAP - NOAA/NCEP Global Forecast System (GFS) Atmospheric Model - Data Access Form

 

여기서 '기온'으로 통상적으로 사용하는 tmp2m(지표 2m 레벨 기온)을 사용했다. 데이터 다운로드 링크를 만들 수 있으므로 웹페이지상에서 최신의 데이터를 실시간으로 fetch 하는 것도 가능하지만, 시간 흐름에 따른 관찰이 어렵고 전송용량이 다소 커서 데이터를 미리 받아서 가공해두기로 했다.

 

데이터는 .nc 포맷으로 다운받았다. 쓸데없이 행 이름등을 반복하지 않고 자체 압축도 지원해서 용량이 매우 작다. 

원본은 위경도 0.5도 간격으로 0~359.5도, 그리고 -90~90도의 그리드 데이터다. 720x361의 해상도를 지닌다.

 

시각화에 사용한 deck.gl의 요구조건에 맞추기 위해 경도를 -180~180도로 조정하고 좌표계 변환을 통해 해상도를 1024x1024로 변환시켰다. 1024x1024는 구글 타일맵 기준 줌레벨 5단계에 해당한다. 전 세계를 EPSG3857 기준으로 가로 32장, 세로32장으로 구분하게 된다. 4단계는 512x512라서 원본대비 경도 데이터의 손실이 있다고 판단했다. 지도를 확대하면 선형 보간을 하면 되기 때문에 굳이 원본의 해상도를 많이 넘어서 6단계 타일맵으로 준비할 필요는 없다.

 

 

PNG 포맷으로 저장

 

래스터 타일맵을 사용하기로 했으므로 PNG 포맷으로 데이터를 저장했는데, 여러단계의 시행착오를 통해 RG, BA채널에 각각 현재 시점과 한단계 뒤의 시점인 3시간 후의 기온값을 담았다. 시행착오 과정은 다음과 같다.

 

1. 첫 번째 시도

RG, BA에 현재시각과 다음 시각을 담음

RG 채널을 예로 들자면, 

R채널 8비트에 기온의 정수값에 100을 더해서 저장. 기온이 156도를 넘어가는 일은 지구 멸망 전에는 없으므로.

G채널 8비트에 기온 실수값을 0~256 해상도로 구분하여 저장. 0.005도보다 조금 더 자세한 정밀도를 지니므로 문제없다고 판단함

 

아래와 같음 그림이 나옴

 

웹에서 PNG를 읽을 때 alpha 값에 따라서 RGB값을 미세하게 보정함. 아예 분명한 방법으로 보정하면 역계산으로  값을 되돌릴 수 있는데 정말 규칙을 종잡을 수가 없음. 테스트해보니 A값이 0이면 모든 RGB값이 0으로 바뀜. A값이 작으면 변동이 심해지는 것으로 판단하여 일단 다른 방법을 모색해보기로 함.

 

2. 두 번째 시도

RGBA 중 A는 못 쓰는 채널이라 생각해서 버리고, 남아있는 8 x 3 = 24 비트를 12비트씩 잘라서 현재시각값과 다음시각값을 넣음

12비트면 4096 단계라서 약간 아쉬워도 그럭저럭 시각화 하는데는 문제가 없다고 판단함

그래서 그렸는데, 현재시각은 문제가 없었는데(없어 보였는데) 다음 시각에서 아래와 같은 그림이 나옴

 

확대하면 아래와 같음

 

원인을 찾느라 한참 걸렸는데, 결과적으로 보간의 보자만 알면 이해할 수 있는 부분을 생각지 못한 어처구니없는 실수였음.

 

준비해놓은 RGB 텍스쳐를 셰이더에서 읽어서 보간했는데, 셰이더는 이미지를 확대할 경우 원본 텍스쳐 픽셀 각각의 채널에 대해 하드웨어적으로 보간함. 그런데 G 채널의 경우 현재시각의 소수점 부분 일부가 상위비트에, 다음시각의 정수 부분 일부가 하위비트에 섞인 채널이므로 이 채널이 보간되는 순간 값은 두 시점의 다른 스케일 데이터가 섞여서 값은 그냥 엉터리가 됨. 위 그림이 그 결과임

 

 

3. 세번째 시도

 

PNG가 아닌 다른 포맷들도 검토해보았으나, pmtiles 자체에서 지원하는 포맷이 한정적이라 결국 PNG에서 해결해야 함.

다시 첫 번째 RG, BA로 나눠 쓰는 방식으로 가기로 함. 대신 알파값에 따른 변형을 최소화하기 위해 알파값에 기본으로 128을 깔고, 나머지 128~255에서만 변위를 주기로 함.

 

그래서 이번에는 첫번째 그림처럼 곰보같은 이미지도 나오지 않고 잘 되는 것처럼 보였으나.....

현재시점 데이터의 다음 시각과, 다음시점 데이터의 현재시각이 일치해야 하지만, 묘하게 1씩 차이나는 부분들이 생김을 발견. 원인은 아까와 같았음. A값이 충분히 크면 RGB가 변형되지 않는다고 생각했는데, 그렇지 않았음. A값이 거의 255라도 변형되는 RGB가 생김.

그런데 정수 부분에 해당하는 R과 B 채널이 1만큼만 바뀌어도 기온이 1도나 올라가버리므로 시간흐름에 따라 자연스럽지 못한 이미지를 만들어내는 것이었음.

 

어찌되었든 이제는 여기서 무조건 해결해야 한다는 생각에 이것저것 상세히 디벼파기 시작함.

같은 라이브러리로 node.js에서 읽으면 제대로 읽는데(우측) 브라우저에서 읽으면 좌측처럼 1씩 변하는 것들이 생김

로딩 옵션에 premultiplyAlpha : none 와 같은 것이 있었는데(로딩시 알파 채널을 미리 곱하지 않는 옵션) 어디까지나 가이드라인같은 옵션이라 브라우저에 따라 지원 안할 수도 있다고 함. 

 

그래서 한참 방법들을 찾아보다가 아래와 같이 디폴트 로더 대신 외부 라이브러리로 png를 읽을 수 있다는 것을 알게 됨

import { decode } from "fast-png";

const loader = ImageLoader;
//@ts-ignore
loader.parse = async (arrayBuffer: ArrayBuffer, options: any) => {
  try {
    const image = await decode(arrayBuffer); // decode가 Promise 반환
    //console.log("Decoded image:", options, image);
    return new ImageData(
      new Uint8ClampedArray(image.data),
      image.width,
      image.height
    );
  } catch (error) {
    console.error("Error decoding image:", error);
    throw error; 
  }
};

 

그 결과 데이터는 왜곡 보정 없이 읽혔음.

결국 첫번째 시도에서 fast-png 라이브러리를 사용하면 되었으나...., 대부분의 시행착오 경과는 이런 식임.

 

데이터를 인코딩하는 코드 부분은 아래와 같음

# R 채널: 정수 부분 (온도에 +100을 더하고 0-255로 클램프)
current_temp = current_temperature_reprojected.data
next_temp = next_temp_reprojected.data

#알파채널이 0으로 가면 손실이 많아져서, 기본적으로 128은 깔고 간다.
current_temp_int = np.clip(np.floor(current_temp + 100), 0, 255).astype(np.uint8)
current_temp_frac = np.clip(np.round(((current_temp + 100) % 1) * 125 + 128), 0, 255).astype(np.uint8)

next_temp_int = np.clip(np.floor(next_temp + 100), 0, 255).astype(np.uint8)
next_temp_frac = np.clip(np.round(((next_temp + 100) % 1) * 125 + 128), 0, 255).astype(np.uint8)

rgba_data = np.stack((current_temp_int, current_temp_frac, next_temp_int, next_temp_frac), axis=-1)


# PIL을 사용하여 이미지 저장
output_file = os.path.join(write_folder, f'{unix_time}-{zoom_level}.png') 
img = Image.fromarray(rgba_data, mode='RGBA')
img.save(output_file, format='PNG')

 

그래서 결과물 중 줌레벨 5단계의 png는 아래와 같음

 

사람의 눈으로 쉽게 이해할 수 있는 이미지는 아니지만, 그래도 육지와 바다의 온도 패턴이 달라서 매직아이처럼 잘 보면 대륙과 바다가 구분되기는 한다.

 

 

 

 

경계를 겹쳐가면서 타일로 만들기

 

 

사실 모든 시행착오는 시간의 전후관계가 얽혀 있지만, 편의상 그냥 하나씩 풀어가는 중이다.

xyz 타일을 만들어 사용했는데, 0~5단계로 미리 준비했다.

0단계는 세계 전체가 32x32 pixel 한장에 담긴다.

1단계는 가로 2등분, 세로 2등분해서 4장에 전 세계가 담긴다.

5단계는 가로 32등분, 세로 32등분해서 1024장에 전 세계가 담긴다.

한 장의 가로세로가 32pixel로 작은 이유는, 우선 특정 줌 단계에서 들어오는 데이터 용량을 최소화 시키려고 했고, 5단계 줌 레벨에서 0.5도 간격의 원본 데이터를 뿌리기로 했다. 그래서 계산해보니 5단계에서 1024x1024의 해상도면 충분했고, 그걸 가로세로 32등분으로 나누니 한 픽셀 기준이 32x32가 되었다.

 

그래서 타일링을 해서 그리고 검증을 위해 1도 간격의 등온선 같은 것을 그려보았는데....

 

자세히 보면 위아래 그림 가운데  타일 경계지점이 묘하게 어긋나 있는 것을 볼 수 있다. 원래 크기인 32x32 타일은 실제로 화면에서는 200x200 이상으로 보일 수도 있는데, 이렇게  타일이 확대되면서 원본 값의 사이사이는 자동적으로 선형 보간을 한다. 그런데 가장자리 같은 경우에는 옆 타일의 값을 모르니 보간하지 못하고 자기 자신을 연장시켜서 가장자리까지 끌고 간다.

 

결국 타일을 만들 때 옆 타일 부분까지 한 줄씩 더 붙여서 만드는 방법밖에는 없다. 그리고는 webgl의 클리핑 기능을 이용해서 잘라서 그리면 옆 타일과 딱! 정확하게 맞아떨어진다.

 

즉, 아래처럼 한 타일을 34x34로 만들어줘야 하고, 가장자리는 34x33 이나 33x34가 될 수 있고, 줌레벨 0인 경우는 32x32가 되어야 한다. 결국 타일이 자연스럽게 그려지도록 하려면 모든 부분들을 커스터마이징해줘야 한다.

 

 

 

그리고 그 '딱 맞아떨어지게' 그리는 부분도, 매우 어려운 난이도는 아니지만, 좌표계 변동이 관련되어 있어 약간 정신차리고 계산을 해야 된다. 이미지는 3857 좌표계라 위아래로 갈수록 위도 1도 차이가 점점 벌어지는데, 클리핑 좌표는 4326으로 넣어야 하기 때문에 좌표 변환을 해가면서 변환해야 한다. 즉, 픽셀 한 줄이 얼마만큼의 거리에 해당하는지 계산해서 타일마다 클리핑 되는 정도를 다르게 지정해야 한다.

그리고 최초 데이터가 0~359.5도까지 존재했으므로 그와 관련된 약간의 offset도 적용해줘야 된다.

function calculateTileOffsets_3857(
  boundingBox: number[][],
  index: { x: number; y: number; z: number },
  tileSize: number,
  overlapPixels: number
) {
  const { x, y, z } = index;

  const R = 6378137; // Earth's radius in meters
  const originShift = (2 * Math.PI * R) / 2; // Half the world's circumference in meters

  // 줌 레벨에서 타일 개수 계산
  const numTiles = Math.pow(2, z);

  // 타일 하나의 크기 (미터 단위)
  const tileExtent = (2 * originShift) / numTiles;

  // 한 픽셀의 크기 계산
  const pixelSize = tileExtent / tileSize;

  // 타일의 경계 좌표 계산 (기본 좌표, 오프셋 전)
  const xMin = -originShift + x * tileExtent;
  const xMax = xMin + tileExtent;
  const yMax = originShift - y * tileExtent;
  const yMin = yMax - tileExtent;

  // 오프셋 계산 (픽셀 크기 * 겹칠 픽셀 수)
  const offset = pixelSize * overlapPixels;

  const leftBottom = proj4("EPSG:3857", "EPSG:4326", [xMin, yMin]);
  const rightTop = proj4("EPSG:3857", "EPSG:4326", [xMax, yMax]);

  const leftBottomOffset = proj4("EPSG:3857", "EPSG:4326", [
    xMin - offset,
    y === numTiles - 1 ? yMin : yMin - offset,
  ]);
  const rightTopOffset = proj4("EPSG:3857", "EPSG:4326", [
    xMax + offset,
    y === 0 ? yMax : yMax + offset,
  ]);


  //최초의 데이터가 0~359.5도까지 있었는데,
  //그 픽셀이 0~360도에 대응하므로, 0.25도만큼 빼준다.
  let adjustOffset = -0.25;

  const boundingBoxRe = [
    leftBottomOffset[0] + adjustOffset,
    leftBottomOffset[1],
    rightTopOffset[0] + adjustOffset,
    rightTopOffset[1],
  ] as [number, number, number, number];

  const clipBounds = [
    leftBottom[0] + adjustOffset,
    leftBottom[1],
    rightTop[0] + adjustOffset,
    rightTop[1],
  ] as [number, number, number, number];

  const boundingBoxWidth = Math.abs(boundingBoxRe[2] - boundingBoxRe[0]);
  // 타일이 각각 줌 레벨의 최상단이거나 최하단에 위치할 경우 상하오프셋을 취소한다.
  if (
    (clipBounds[2] - adjustOffset) * (boundingBoxRe[2] - adjustOffset) < 0 &&
    boundingBoxWidth > 90
  )
    boundingBoxRe[2] += 360;
  if (
    (clipBounds[0] - adjustOffset) * (boundingBoxRe[0] - adjustOffset) < 0 &&
    boundingBoxWidth > 90
  )
    boundingBoxRe[0] -= 360;

  return {
    boundingBoxRe,
    clipBounds,
  };

 

 

 

겹치는 부분을 확인해보기 위해 클리핑을 하지 않고 반투명하게 그리면 타일이 겹치는 부분들이 아래처럼 보인다.

잘 보면 각각의 경계 부분 값이 가장자리까지 연장되어 있고, 겹치는 부분이 정확하게 일치하는 것을 확인할 수 있다.

 

다행히 deck.gl에서 ClipExtension으로 클리핑 기능을 기본 지원해주는 덕분에, 클리핑 하는 과정은 원래의 파이프라인을 건드리지 않고 해결되었다.

 

 

타일 데이터 관리

 

결론부터 말해보자면, 일단 ERDDAP가 공개하는 2022년 12월 1일부터 현재 시점까지 2년 x 365일 x 3시간 간격 8개 시점 = 약 5800개의 Pmtiles 형식의 타일 데이터를 만들었고, 하나하나는 3.5MB정도 된다. 이 파일들을 Github의 공개 저장소에 올렸다. 공개 저장소에 올릴 경우 다른 사이트에서도 받아갈 수 있다. CORS 제한없이 허용된다. 

github의 경우 5GB 이상은 권장되지 않는다고 하는데, 혹시라도 관련된 이메일을 깃헙 측으로부터 받으면 그 때가서 고민해보기로 했다. 저장소 20GB가 요새의 관점에서 그렇게 큰 용량은 아니니.

 

한 시점 당 1개의 Pmtiles 를 만든 것인데, 만약 그냥 폴더 구조의 xyz 타일을 만들 경우 5800개 시점 x 1365 = 7,917,000 개의 파일을 관리해야 한다. 복사부터 액세스 속도까지 이것저것 모두 문제가 된다. 그런 면에서 Pmtiles는 최초에 인덱스 헤더 16kb 정도를 받아야 하는 단점은 있지만, 그 다음부터는 3ms 정도의 매우 빠른 응답속도를 보여준다. 

처음에 폴더구조의 xyz으로 테스트해봤을 때, 15ms 정도의 느린 응답을 보여주었었다. 즉 pmtiles가 4배 정도 빠르다.

 

단, 처음에 pmtiles 400여개 올렸을 때는 인덱스도 20ms 를 채 넘지 않는 시간에 넘어왔던 것 같은데, 5800개로 많아지니 인덱스 부분 받아오는데 1초 내외가 걸린다. 5800개를 100여개 폴더에 분할해서 넣어도 마찬가지다. localhost로 ssd 에서 테스트 할 때는 당연히 인덱스든 타일이든 3ms 정도에 모두 받아온다.

유료 저장소로 바꾸면 좀 나아질지도 모르겠는데, 그런건 나중에 정말 필요할 때 찾아서 테스트 해보기로 했다.

 

 

 

여기까지가 데이터에 대한 얘기다. 이제 그리는 이야기를 약간 해보자.

 

 

 

타일 데이터로 전송 패킷 경량화

 

 

타일 데이터를 사용한 이유는 위의 GIF에서 확인할 수 있다. 타일 경계마다 하이라이트 되도록 했다.

 

특정 시점의 한 화면을 그리는데 4장~9장 정도의 타일만 필요하다. 즉 30KB 내외로 한 시점을 그릴 수 있다는 것. 

속도 측면에서 볼 때 파일 경량화는 끊김없는 애니메이션 뷰를 만들어주는데 필수적이다. 그리고 지구 온난화 시대에 기온 데이터를 시각화 하는데 웹 트래픽이나 저장공간을 되도록 줄여야 하는게 아닐까 하는 생각도 했다. (사실 지구 온난화 시대에 진정으로 기여하려면 이런 시각화 자체를 아예 만들지 않아서 전기와 자원 사용량을 줄이는게 나을 것 같기도 하다)

 

단점도 있다. 1단계까지 축소되면 전세계가 64x64 의 데이터로 표현되는데, 그만큼 해상도가 낮은 데이터 기반으로 보간을 하게 되므로 아래 그림처럼 축소 상태에서는 측정 지점의 기온이 변할 수 밖에 없다. 물론 넓은 시점에서의 기온 패턴이 변하는 것은 아니지만 국부적으로는 기온값이 틀린것처럼 보일 수 있다는 말이다.

 

 

 

@loader.gl/pmtiles

 

아직 deck.gl에 통합되지는 않았지만, loaders.gl에는 pmtiles를 읽어올 수 있는 모듈이 있다.

해당 모듈로 TileLayer에서 타일 데이터를 읽어서  deck.gl의 서브레이어에 전달해주었다.

 

파일마다 헤더에 들어있는 인덱스를 기반로 동작하는 pmtiles의 구조를 약간 알아야 구현할 수 있지만, 다른 작업을 하느라 golang 원본을 Java로 포팅해본 경험이 있어, 다행히 어렵지 않게 구현하였다.

 

 

 

TileLayer와 subrender로서의 CustomBitmapLayer

 

최대한 기존 라이브러리를 살려보자는 생각에, Deck.gl의 TileLayer를 사용해서 작업했다.

기본적인 TileLayer의 동작은, Tile로 받아온 데이터를 subrender 메서드를 통해 bitmaplayer 등으로 표현하는 방식이다. BitmapLayer를 상속해서 CustomBitmapLayer를 만들어 사용했다.

 

프래그먼트 레이어는 거의 새로 쓰다시피 했으며, uniform 변수도 몇개 더 넣어주었다. 태양의 고도를 계산하려면 unixtime 등 추가적인 변수가 셰이더로 전달되어야 한다.

 

 

표출 온도 필터링

 

셰이더의 코드는 GPU가 곧바로 실행하기 때문에 어떤 작업들은 굉장히 빠르게 반응한다.

기본적으로 온도에 따라 그라데이션을 주거나, smoothStep을 이용해서 5도 간격의 등온선을 그리거나, 태양의 고도에 따라 어두운 색상을 오버레이 하는 작업들이 모두 20-30줄의 코드 안에서 수행된다.

 

그리고 uniform 변수를 통해 최대최소 표출 온도 한계값을 입력하면 아래 그림처럼 보고자 하는 온도 범위만 볼 수 있다.

 

데이터 시각화에서는 이런 간단한 필터링 장치들이 어떤 현상을 관찰하는데 매우 도움이 된다.

 

 

아래 그림에서는 태양 고도 60도 이상인 부분(가장 밝게 보이는 부분)과 기온이 상승하는 지역의 상관성을 직관적으로 관찰할 수 있다.

 

 

 

 

시간 변화에 따른 자연스러운 기온 변화

 

위와 같은 기나긴 시행착오를 겪으면서 그려보려고 했던 것은 아래와 같은 '자연스러운' 애니메이션이다.

 

windy.com 이나 zoom.earth, earth.nullschool.net 모두 많은 데이터를 볼 수는 있는데 시간 변화를 자연스럽게 보여주지는 못한다. 물론 windy.com이나 zoom.earth의 경우 낮은 줌레벨에서도 고해상도의 데이터를 불러오느라 그랬을 것 같기는 한데(게다가 애니메이션 용으로 데이터를 방출해주면 클라우드 비용도 만만치 않을테니...), 여튼 애니메이션을 만들기 위해 극복해야 할 여러가지 기술적인 목표들을 설정하고 하나씩 해결해 나간 부분이 의미있는 일이었다. 

 

 

 

 

나가며

 

사람이 느끼는 쾌적한 실외 온도인 15~25도로 필터링해놓고 세계지도를 보면, 대부분의 계절에서 그 혜택을 길게 받는 지역이 별로 없다. 그나마 우리나라는 아슬아슬한 경계에 걸쳐져 있어서 봄과 가을이 존재했는데, 사실 올해는 가을다운 가을이 거의 사라졌었던 것 같다. 여름에서 갑자기 겨울이 된 기분이다.

 

이 시각화로 기후위기를 곧바로 알 수는 없지만, 시간의 전후관계 속에서 관찰을 하면 우리나라 기온이 주변 수온과 북쪽의 찬 기운으로부터 영향받는다는 것을 '어렴풋이' 볼 수 있다. 아마도 대기과학의 전문가라면 이런 즉각적인 반응성을 지닌 툴을 통해서 설명해 줄 수 있는 부분이 좀 더 많을 것이라 생각한다.

대기과학자들과 작업할 일이 있을지는 모르겠지만, 그럴 일이 있다면 '어떠한 시각적 장치를 더 넣어서 안 보이던 것을 잘 보이게 해 줄 것인가'에 대해 좀 더 생각할 기회가 되겠지. 

 

 

"데이터 시각화가 무엇을 하는 일입니까"라고 누군가 물어올 때, 나는 요새 "보이지 않던 것을 보이게 해 줍니다"라고 대답한다.

현미경이 안 보이던 세상을 확대해서 보이게 해 주는 것이라고 한다면, 이번에 만든 시각화는 멈춰 있던 세상을 앞뒤로 움직이게 해주는 것이라고나 할까.

 

 

 

기후 데이터 뷰어