junyeokk
Blog
Parallel Computing·2025. 10. 23

3D 배열과 이미지 필터링

3D 데이터 처리의 필요성

GPU가 3D 데이터를 처리하는 이유

300 × 300 × 256 크기의 3D 이미지 필터링을 구현했을 때 놀라운 결과가 나타났다. CPU에서 44,784 usec가 걸렸지만 GPU에서는 830 usec로 54배 빠르게 처리했다. 같은 연산인데 왜 GPU가 이렇게 빠를까?

답은 데이터의 구조에 있다. 3D 배열은 공간의 모든 점에서 독립적인 계산을 수행한다. 300 × 300 × 256 = 23,040,000개의 원소가 있고, 각 원소는 다른 원소와 무관하게 필터링된다. CPU는 이 계산을 순차적으로 하나씩 수행하지만, GPU는 수천 개의 스레드로 동시에 처리한다. 병렬 처리가 가능한 데이터 구조가 GPU의 힘을 극대화한다.

그렇다면 3D 데이터는 어디에 사용될까? 의료 데이터라면 256개의 2D 평면이 쌓인 구조다.

하지만 실제 메모리는 1D로 배치된다. 컴퓨터 메모리는 주소 0부터 시작하는 일직선 공간이므로, 3D 좌표를 1D 주소로 변환해야 한다.

plain
i = (z * HEIGHT + y) * WIDTH + x;

이 공식은 row-major 방식의 인덱스 계산이다. z=0 평면 전체가 먼저 배치되고, z=1 평면이 그 뒤를 따른다. 각 평면 내에서는 2D 행렬처럼 행 단위로 배치된다. 예를 들어 M[0][0][0], M[0][0][1], M[0][0][2], ..., M[0][1][0], M[0][1][1], ... 순서로 메모리에 저장된다.

CT 스캔 데이터가 이렇게 저장된다. 환자의 머리부터 발끝까지 촬영한 2D 이미지를 순서대로 메모리에 쌓으면 3D 볼륨이 완성된다. X-ray 1장이 2D 행렬 하나가 되고, 256장을 쌓으면 3D 행렬이 된다.

CT 스캐닝이란?
Computed Tomography의 약자로 컴퓨터 단층 촬영을 의미한다. X-ray를 여러 각도에서 촬영하여 신체의 단면 이미지를 생성하고, 이를 쌓아서 3D 구조를 재구성하는 의료 영상 기술이다.

이미지 필터링의 원리

이미지 필터링은 원본 이미지의 각 픽셀에 필터 값을 곱하여 새로운 이미지를 만드는 연산이다. element-by-element multiplication이라고도 부른다.

plain
C[i] = A[i] * B[i]

A는 원본 이미지, B는 필터, C는 결과 이미지다. 각 픽셀끼리 곱한다. 가우시안 필터를 적용하면 흐릿한 효과가 나타나고, 샤프닝 필터를 적용하면 경계가 강조된다. 3D 의료 영상에서는 노이즈를 제거하거나 특정 조직을 강조하는 데 사용한다.

기본 3D 필터링 구현

CPU에서 3D 필터링 구현하기

300 × 300 × 256 크기의 3D 이미지 필터링을 CPU에서 구현해보자. 23,040,000개의 원소를 하나씩 처리한다.

c
dim3 dimImage(300, 300, 256); // width, height, depth

int main(const int argc, const char* argv[]) {
    matA = new float[dimImage.z * dimImage.y * dimImage.x];
    matB = new float[dimImage.z * dimImage.y * dimImage.x];
    matC = new float[dimImage.z * dimImage.y * dimImage.x];

    ELAPSED_TIME_BEGIN(0);
    for (register unsigned z = 0; z < dimImage.z; ++z) {
        for (register unsigned y = 0; y < dimImage.y; ++y) {
            for (register unsigned x = 0; x < dimImage.x; ++x) {
                int i = (z * dimImage.y + y) * dimImage.x + x;
                matC[i] = matA[i] * matB[i];
            }
        }
    }
    ELAPSED_TIME_END(0);
}

위 코드는 3중 루프로 모든 원소를 순회한다. 가장 안쪽 루프가 x축(너비)이므로 메모리를 순차적으로 접근한다. 각 루프에서 (z * dimImage.y + y) * dimImage.x + x 공식으로 3D 인덱스를 1D 인덱스로 변환한다. z=0 평면 전체를 먼저 처리하고, z=1 평면으로 넘어간다.

Intel Core i5-3570에서 44,784 usec가 소요된다. 약 45밀리초다.

GPU에서 3D 필터링 구현하기

GPU에서는 3D thread block으로 3D 데이터를 처리한다. 8 × 8 × 8 = 512개 스레드로 구성된 블록을 사용한다. 300 × 300 × 256 크기를 커버하려면 300 / 8 = 37.5이므로 38개 블록이 필요하다. 256 / 8 = 32개 블록이 필요하다. grid 크기는 38 × 38 × 32가 된다.

c
dim3 dimBlock(8, 8, 8);
dim3 dimGrid(div_up(dimImage.x, dimBlock.x),
             div_up(dimImage.y, dimBlock.y),
             div_up(dimImage.z, dimBlock.z));

div_up 함수는 올림 나눗셈을 수행한다. (lhs + rhs - 1) / rhs로 구현한다. 300을 8로 나누면 37.5인데 정수 나눗셈으로는 37이 나온다. div_up을 사용하면 38이 나온다. 나머지가 있으면 자동으로 올림된다.

왜 (lhs + rhs - 1) / rhs로 올림을 구현할까?
정수 나눗셈에서 올림을 구현하는 수학적 트릭이다. lhs = 10, rhs = 3이면 (10 + 3 - 1) / 3 = 12 / 3 = 4가 된다. 일반 나눗셈 10 / 3 = 3보다 1 크다. lhs가 rhs의 배수라면 올림이 일어나지 않고, 나머지가 있으면 자동으로 1이 추가된다.

c
template <typename TYPE>
__host__ __device__ inline TYPE div_up(TYPE lhs, TYPE rhs) {
    return (lhs + rhs - 1) / rhs;
}

커널 함수는 각 스레드가 하나의 원소를 담당한다.

c
__global__ void kernel_filter(float* c, const float* a, const float* b,
                               unsigned ndim_z, unsigned ndim_y, unsigned ndim_x) {
    unsigned idx_z = blockIdx.z * blockDim.z + threadIdx.z;
    unsigned idx_y = blockIdx.y * blockDim.y + threadIdx.y;
    unsigned idx_x = blockIdx.x * blockDim.x + threadIdx.x;

    if (idx_x < ndim_x && idx_y < ndim_y && idx_z < ndim_z) {
        unsigned i = (idx_z * ndim_y + idx_y) * ndim_x + idx_x;
        c[i] = a[i] * b[i];
    }
}

각 스레드는 blockIdx와 threadIdx를 사용해서 자신이 처리할 3D 좌표 (idx_x, idx_y, idx_z)를 계산한다. 범위를 벗어나는 스레드는 if 조건으로 걸러진다. grid가 38 × 8 = 304인데 실제 크기는 300이므로 4개 스레드는 아무 일도 하지 않는다. 3D 좌표를 1D 인덱스로 변환하여 배열에 접근한다.

커널 호출은 다음과 같다.

c
kernel_filter <<<dimGrid, dimBlock>>> (dev_matC, dev_matA, dev_matB,
                                        dimImage.z, dimImage.y, dimImage.x);

GeForce RTX 2070에서 830 usec가 소요된다. CPU 대비 54배 빠른 속도다.

3D Pitched Memory 사용

3D Pitched Memory란 무엇인가

2D 배열에서 cudaMallocPitch를 사용했던 것처럼 3D 배열에서도 pitched memory를 사용할 수 있다. 일반 3D 배열은 WIDTH를 사용하지만, pitched 3D 배열은 PITCH를 사용한다.

plain
일반:   index = (z * HEIGHT + y) * WIDTH + x
pitched: index = (z * HEIGHT + y) * PITCH + x

WIDTH는 실제 데이터 너비지만, PITCH는 메모리 정렬을 위해 확장된 너비다. 300 × 300 × 256 크기에서 WIDTH = 300 × sizeof(float) = 1200바이트이지만, PITCH는 512바이트 경계로 올림되어 1536바이트가 될 수 있다. 336바이트가 낭비되지만 메모리 접근이 정렬되어 성능이 향상된다.

정렬 경계는 시스템에 따라 다르다. 32B, 64B, 128B, 256B, 512B, 1024B 등 다양하다. CUDA 버전과 하드웨어에 따라 변한다. 프로그래머가 이 값을 추측할 필요는 없다. cudaMalloc3D 함수가 자동으로 계산해서 알려준다.

cudaMalloc3D와 cudaPitchedPtr

cudaMalloc3D는 3D pitched memory를 할당한다.

c
__host__ cudaError_t cudaMalloc3D(cudaPitchedPtr* pitchedDevPtr, cudaExtent extent);

extent로 원하는 크기를 지정하면, 함수가 메모리를 할당하고 pitchedDevPtr에 결과를 저장한다. cudaPitchedPtr 구조체는 할당된 메모리 정보를 담는다.

c
struct cudaPitchedPtr {
    void*  ptr;    // 실제 메모리 주소
    size_t pitch;  // 바이트 단위 pitch
    size_t xsize;  // 바이트 단위 너비
    size_t ysize;  // 원소 단위 높이
};

300 × 300 × 256 배열을 할당한다고 하자. xsize = 300 × sizeof(float) = 1200바이트다. pitch는 시스템이 계산해서 1536바이트로 설정할 수 있다. 512바이트 경계로 올림한 값이다. ysize = 300이다.

호스트 메모리는 일반 malloc으로 할당했으므로 pitch가 없다. make_cudaPitchedPtr로 pitched 포인터 형태로 만들 수 있다.

c
struct cudaPitchedPtr pitchedA
    = make_cudaPitchedPtr(matA, dimImage.x * sizeof(float),
                          dimImage.x * sizeof(float), dimImage.y);

호스트 메모리는 연속 배치되므로 pitch와 xsize가 같다.

cudaMemcpy3D로 복사하기

3D pitched memory 간의 복사는 cudaMemcpy3D를 사용한다. cudaMemcpy3DParams 구조체로 복사 정보를 전달한다.

c
struct cudaMemcpy3DParams {
    struct cudaPos srcPos;          // 소스 시작 위치
    struct cudaPitchedPtr srcPtr;   // 소스 pitched 포인터
    struct cudaPos dstPos;          // 목적지 시작 위치
    struct cudaPitchedPtr dstPtr;   // 목적지 pitched 포인터
    struct cudaExtent extent;       // 복사할 크기
    enum cudaMemcpyKind kind;       // 복사 방향
};

호스트 메모리를 준비하고 pitched 포인터를 만든다.

c
matA = new float[dimImage.z * dimImage.y * dimImage.x];
struct cudaPitchedPtr pitchedA
    = make_cudaPitchedPtr(matA, dimImage.x * sizeof(float),
                          dimImage.x * sizeof(float), dimImage.y);

디바이스 메모리를 할당한다. make_cudaExtent로 크기를 지정한다.

c
struct cudaExtent extentInByte = make_cudaExtent(dimImage.x * sizeof(float),
                                                  dimImage.y, dimImage.z);
struct cudaPitchedPtr dev_pitchedA = { 0 };
cudaMalloc3D(&dev_pitchedA, extentInByte);

cudaExtent는 3D 크기를 나타낸다. make_cudaExtent(w, h, d)에서 w는 바이트 단위 너비, h는 원소 단위 높이, d는 원소 단위 깊이다. 300 × 300 × 256 배열이면 make_cudaExtent(300 × sizeof(float), 300, 256)이 된다.

복사 매개변수를 설정하고 cudaMemcpy3D를 호출한다.

c
struct cudaPos pos_origin = make_cudaPos(0, 0, 0);
struct cudaMemcpy3DParams paramA = { 0 };
paramA.srcPos = pos_origin;
paramA.srcPtr = pitchedA;
paramA.dstPos = pos_origin;
paramA.dstPtr = dev_pitchedA;
paramA.extent = extentInByte;
paramA.kind = cudaMemcpyHostToDevice;
cudaMemcpy3D(&paramA);

srcPos와 dstPos는 복사 시작 위치다. 전체 배열을 복사하므로 (0, 0, 0)이다. kind는 cudaMemcpyHostToDevice로 호스트에서 디바이스로 복사한다. cudaMemcpy3D는 각 평면을 개별적으로 복사하면서 pitch 차이를 자동으로 처리한다.

3D Pitched 커널 구현하기

pitched memory를 사용하는 커널은 pitch 매개변수를 받는다. 인덱스 계산에서 WIDTH 대신 pitch를 사용한다.

c
__global__ void kernel_filter(void* matC, const void* matA, const void* matB,
                               size_t pitch, unsigned ndim_z, unsigned ndim_y, unsigned ndim_x) {
    register unsigned idx_z = blockIdx.z * blockDim.z + threadIdx.z;
    register unsigned idx_y = blockIdx.y * blockDim.y + threadIdx.y;
    register unsigned idx_x = blockIdx.x * blockDim.x + threadIdx.x;

    if (idx_x < ndim_x && idx_y < ndim_y && idx_z < ndim_z) {
        register unsigned offset_in_byte = (idx_z * ndim_y + idx_y) * pitch + idx_x * sizeof(float);
        *((float*)((char*)matC + offset_in_byte))
            = *((const float*)((const char*)matA + offset_in_byte))
            * *((const float*)((const char*)matB + offset_in_byte));
    }
}

offset_in_byte는 바이트 단위로 계산한다. (idx_z * ndim_y + idx_y) * pitch는 해당 행의 시작 위치(바이트)를 가리킨다. idx_x * sizeof(float)는 그 행에서 idx_x 번째 원소의 오프셋(바이트)이다. char* 캐스팅으로 바이트 단위 연산을 수행한 뒤 float*로 다시 캐스팅한다.

커널 호출 시 dev_pitchedA.ptr로 실제 포인터를, dev_pitchedA.pitch로 pitch 값을 전달한다.

c
kernel_filter <<<dimGrid, dimBlock>>>(dev_pitchedC.ptr, dev_pitchedA.ptr, dev_pitchedB.ptr,
                                       dev_pitchedA.pitch, dimImage.z, dimImage.y, dimImage.x);

GeForce RTX 2070에서 844 usec가 소요된다.

성능 비교

300 × 300 × 256 3D 이미지 필터링 성능을 정리하면 다음과 같다.

plain
CPU:                     44,784 usec
GPU (일반):                 830 usec (54배 빠름)
GPU (pitched):              844 usec (53배 빠름)

GPU는 CPU 대비 54배 빠르다. pitched 버전은 일반 버전과 거의 비슷한 성능을 보인다. 이 테스트에서는 pitched memory의 장점이 크게 드러나지 않았지만, 더 복잡한 메모리 접근 패턴에서는 차이가 커진다. 메모리 정렬은 일관된 성능을 보장하는 안전장치다.

관련 문서