Indexação do kernel 3D CUDA para filtragem de imagens?

Eu tenho uma matriz de recurso de imagem A é n * m * 31 matriz preenchida para filtragem e tenho B como um filtro de object k * l * 31 . Eu quero obter uma matriz de saída C é p * r * 31 com o tamanho da imagem A sem preenchimento . Eu tento escrever um código CUDA para executar o filtro B sobre A e obter C.

Eu suponho que cada operação de filtragem sobre A com filtro B ocupada por um bloco de threads, assim haverá uma operação k * l dentro de cada bloco de threads. E cada operação de filtragem deslocada será realizada em um bloco de threads diferente. Para A (0,0) a filtragem será em thread_block (0,0) e para A (0,1) será em thread_block (1,0), e assim por diante. Também tenho terceira dimensão como 31. Cada espaço na terceira dimensão será calculado em si. Portanto, com a indexação 3D correta sobre as matrizes, talvez eu possa ter toda a operação de forma muito paralela.

Então a operação é

A n*m*31 XB k*l*31 = C p*r*31 

Como eu poderia fazer a indexação do kernel para as operações de forma eficiente?

Aqui está algum código que escrevi para fazer mais ou menos o que você descreve.

Ele cria um dataset 3D (chamado cell , mas seria comparável ao seu array A ), preenche-o com dados randoms e calcula um array de saída 3D resultante (chamado node , mas seria comparável ao array C ) nos dados em. O tamanho dos dados de A é maior que o tamanho de C (“preenchido” como você o chama) para permitir a passagem da function B sobre os elementos de limite de A No meu caso, a function B é simplesmente encontrar o mínimo dentro do volume cúbico 3D de A que está associado ao tamanho de B (o que chamo de WSIZE , criando uma região 3D WSIZE x WSIZE x WSIZE ) e armazenando o resultado em C

Esse código em particular tenta explorar a reutilização de dados copiando uma determinada região da input A na memory compartilhada de cada bloco. Cada bloco calcula múltiplos pontos de saída (isto é, calcula B um número de vezes para preencher uma região de C) de modo a explorar a oportunidade de reutilização de dados de cálculos vizinhos de B.

Isso pode ajudar você a começar. Você obviamente teria que replace B (meu código de localização mínima) por qualquer que seja a function B desejada. Além disso, você precisaria modificar o domínio B de um cúbico no meu caso para qualquer tipo de prisma retangular correspondente às suas dimensões B. Isso também afetará a operação de memory compartilhada, portanto, convém dispensar a memory compartilhada para sua primeira iteração apenas para obter as coisas funcionalmente corretas e, em seguida, include a otimização de memory compartilhada para ver quais benefícios você pode obter.

 #include  #include  // these are just for timing measurments #include  // Computes minimum in a 3D volume, at each output point // To compile it with nvcc execute: nvcc -O2 -o grid3d grid3d.cu //define the window size (cubic volume) and the data set size #define WSIZE 6 #define DATAXSIZE 100 #define DATAYSIZE 100 #define DATAZSIZE 20 //define the chunk sizes that each threadblock will work on #define BLKXSIZE 8 #define BLKYSIZE 8 #define BLKZSIZE 8 // for cuda error checking #define cudaCheckErrors(msg) \ do { \ cudaError_t __err = cudaGetLastError(); \ if (__err != cudaSuccess) { \ fprintf(stderr, "Fatal error: %s (%s at %s:%d)\n", \ msg, cudaGetErrorString(__err), \ __FILE__, __LINE__); \ fprintf(stderr, "*** FAILED - ABORTING\n"); \ return 1; \ } \ } while (0) // device function to compute 3D volume minimum at each output point __global__ void cmp_win(int knode[][DATAYSIZE][DATAXSIZE], const int kcell[][DATAYSIZE+(WSIZE-1)][DATAXSIZE+(WSIZE-1)]) { __shared__ int smem[(BLKZSIZE + (WSIZE-1))][(BLKYSIZE + (WSIZE-1))][(BLKXSIZE + (WSIZE-1))]; int tempnode, i, j, k; int idx = blockIdx.x*blockDim.x + threadIdx.x; int idy = blockIdx.y*blockDim.y + threadIdx.y; int idz = blockIdx.z*blockDim.z + threadIdx.z; if ((idx < (DATAXSIZE+WSIZE-1)) && (idy < (DATAYSIZE+WSIZE-1)) && (idz < (DATAZSIZE+WSIZE-1))){ smem[threadIdx.z][threadIdx.y][threadIdx.x]=kcell[idz][idy][idx]; if ((threadIdx.z > (BLKZSIZE - WSIZE)) && (idz < DATAZSIZE)) smem[threadIdx.z + (WSIZE-1)][threadIdx.y][threadIdx.x] = kcell[idz + (WSIZE-1)][idy][idx]; if ((threadIdx.y > (BLKYSIZE - WSIZE)) && (idy < DATAYSIZE)) smem[threadIdx.z][threadIdx.y + (WSIZE-1)][threadIdx.x] = kcell[idz][idy+(WSIZE-1)][idx]; if ((threadIdx.x > (BLKXSIZE - WSIZE)) && (idx < DATAXSIZE)) smem[threadIdx.z][threadIdx.y][threadIdx.x + (WSIZE-1)] = kcell[idz][idy][idx+(WSIZE-1)]; if ((threadIdx.z > (BLKZSIZE - WSIZE)) && (threadIdx.y > (BLKYSIZE - WSIZE)) && (idz < DATAZSIZE) && (idy < DATAYSIZE)) smem[threadIdx.z + (WSIZE-1)][threadIdx.y + (WSIZE-1)][threadIdx.x] = kcell[idz+(WSIZE-1)][idy+(WSIZE-1)][idx]; if ((threadIdx.z > (BLKZSIZE - WSIZE)) && (threadIdx.x > (BLKXSIZE - WSIZE)) && (idz < DATAZSIZE) && (idx < DATAXSIZE)) smem[threadIdx.z + (WSIZE-1)][threadIdx.y][threadIdx.x + (WSIZE-1)] = kcell[idz+(WSIZE-1)][idy][idx+(WSIZE-1)]; if ((threadIdx.y > (BLKYSIZE - WSIZE)) && (threadIdx.x > (BLKXSIZE - WSIZE)) && (idy < DATAYSIZE) && (idx < DATAXSIZE)) smem[threadIdx.z][threadIdx.y + (WSIZE-1)][threadIdx.x + (WSIZE-1)] = kcell[idz][idy+(WSIZE-1)][idx+(WSIZE-1)]; if ((threadIdx.z > (BLKZSIZE - WSIZE)) && (threadIdx.y > (BLKYSIZE - WSIZE)) && (threadIdx.x > (BLKXSIZE - WSIZE)) && (idz < DATAZSIZE) && (idy < DATAYSIZE) && (idx < DATAXSIZE)) smem[threadIdx.z+(WSIZE-1)][threadIdx.y+(WSIZE-1)][threadIdx.x+(WSIZE-1)] = kcell[idz+(WSIZE-1)][idy+(WSIZE-1)][idx+(WSIZE-1)]; } __syncthreads(); if ((idx < DATAXSIZE) && (idy < DATAYSIZE) && (idz < DATAZSIZE)){ tempnode = knode[idz][idy][idx]; for (i=0; i>>(d_node, d_cell); cudaCheckErrors("Kernel launch failure"); // copy output data back to host cudaMemcpy(node, d_node, ((nx*ny*nz)*sizeof(int)), cudaMemcpyDeviceToHost); cudaCheckErrors("CUDA memcpy3 failure"); t2 = clock(); t2sum = ((double)(t2-t1))/CLOCKS_PER_SEC; printf(" Device compute took %3.2f seconds. Beginning host compute.\n", t2sum); // now compute the same result on the host for (u=0; u cell[i+u][j+v][k+w]) temphnode = cell[i+u][j+v][k+w]; hnode[u][v][w] = temphnode; } t3 = clock(); t3sum = ((double)(t3-t2))/CLOCKS_PER_SEC; printf(" Host compute took %3.2f seconds. Comparing results.\n", t3sum); // and compare for accuracy for (i=0; i