http://peghoty.blog.163.com/blog/static/493464092013016113254852/
http://blog.csdn.net/augusdi/article/details/12833235
CUDA存储器模型：http://blog.csdn.net/endlch/article/details/44538801
CUDA限定符：http://blog.csdn.net/shouhuxianjian/article/details/42427285
思想即是将内存数据拷贝到显存，在显存上执行并行运算，将结果数据从显存拷贝回内存。
CUDA内有thrust库，类似于C++ stl库。
===========以下是原文=========
挖坑待填。
以上是本机CUDA参数。
需要了解的概念：线程束(wrap)，共享内存，常量内存，纹理内存(?，图形学相关，略)，流，原子操作。
寄存器
寄存器是GPU片上高速缓存， 执行单元可以以极低的延迟访问寄存器。寄存器的基本单元式寄存器文件，每个寄存器文件大小为32bit。局部存储器对于每个线程，局部存储器也是私有的。如果寄存器被消耗完。数据将被存储在局部存储器中。如果每个线程使用了过多的寄存器，或声明了大型结构体或数据，或者编译器无法确定数据的大小，线程的私有数据就有可能被分配到local memory中，一个线程的输入和中间变量将被保存在寄存器或者是局部存储器中。局部存储器中的数据被保存在显存中，而不是片上的寄存器或者缓存中，因此对local memory的访问速度很慢。
共享存储器
共享存储器(share memeory)也是GPU片内缓存存储器。它是一块可以被同一block中的所有线程访问的可读存储器。
使用关键字share添加到变量的声明中，这将使这个变量驻留在共享内存中。cuda c编译器对共享内存中的变量与普通变量将采取不同的处理方式。
对于在GPU上启动的每个线程块，cuda c编译器都将创建该变量的一个副本，线程块中的每一个线程都共享这块内存，但这个线程却无法看到也不能修改其他线程块的变量的副本。这就实现了一种非常好的方式，使得一个线程块中的多个线程能够在计算上进行通信和协作，而且，共享内存缓冲区驻留在物理GPU上，而不是驻留在GPU之外的系统内存中。
常量内存
__constant__将把变量的访问限制为只读。
在接受了这种限制之后，我们希望得到某种回报，与全局内存中读数据相比，从常量内存中读取相同的数据可以节约内存的带宽，主要有两个原因：
-对常量内存的单次读操作可以广播到其他的“领进”线程，这将节约15次读取操作。
-常量内存的数据缓存起来，因此对相同地址的连续读取操作将不会产生额外的内存通信量。
“邻近”是指半个warp中的线程。当处理常量内存时。nvidia硬件将把单次内存读取操作广播到每个半线程束。在半线程束中包含了16个线程，即线程束中数量的一半。
如果在半线程束中的每一个线程访问相同的常量内存地址。那么GPU只会发生一次读操作事件并在随后将数据广播到每个线程。
如果从常量内存中读取大量的数据，那么这种方式产生的内存流量只是全局内存时的1/16.然而，当使用常量内存时也可能产生负面影响。
如果半线程束的所有16个线程需要访问常量内存中不同的数据，那么这个16次读取操作会被串行化，从而需要16倍的时间来发出请求。
但如果从全局内存中读取，那么这些请求会同时发出。在这种情况下，从常量内存读取就慢于从全局内存中读取。
全局存储器
全局存储器(global memeory)位于显存（占据了大部分的显存)。
整个网格中的任意线程都能读写全局存储器的任意位置。在目前的架构中，全局存储器没有缓存。
=======补======
SM与线程束
通常，线程块的数量为GPU中央处理器数量的2倍时，将达到最优性能。
GPU拥有数百个核，其中，SM代表多流处理器，即计算核心，而每个SM又包含8个标准流处理器SP，以及其他。
隶属于同一个SM的8个SP共用同一套取指与发射单元，共用一块共享存储器。
kernel以block为单位执行。
一个block必须被分配到同一块SM中，block的每个thread会发送到一个SP上执行，但一个SM中同一时刻可以有多个活动线程块在等待执行。
同一个block的thread开始于相同的指令地址，理论上能够按不同的分支执行，但实际上由于8个SP共用一套取值与发射单元，同一warp的线程执行的指令是相同的。
如果一个warp的线程跳转如分支语句的同一分支，那么实际执行时间就是这个分支执行时间；
否则，SM需要把每一个分支的指令发射到每个SP，执行时间是执行多个分支的所用时间之和。
故CUDA程序尽量避免分支，尽量warp内不分支。
线程束(warp)：一个线程束由32个连续的线程组成。（简单地说，warp包含32个线程是因为每发射一条warp指令，SM中的8个SP就会将这条指令执行4遍）。warp才是真正的执行单位。
杂
原子操作。同时对global内存写操作，可分批进行，改成先线程块对shared内存写操作，结束后shared内存写入global内存。
__syncthreads()实现了线程块内的线程同步，当任意线程运行到BAR标记处后，暂停运行，直到整个block中所有的thread都运行到BAR标记处后才继续执行。
__syncthreads()勿置于分支语句内。
流：名义上多个流，实际上可能就是kenel（GPU运算）和copy（CPU与GPU间数据传送）两个流。每个流都是一个队列，事件被push入队列中等待执行。
for循环的任务切分的时候，有两种方式划分任务。
1.划分成k段，每段由某个线程执行。
2.按模k同余进行划分，for循环每次的递增量为块大小。
一般第2种方式更优，因为是并行执行，故第二种方式保证每次执行的的时候，不同线程访问的数据位置相邻。
并行算法
归约运算: 每次折半。以求和为例，第一次前1/2 + 后1/2；第二次 前1/4 + 后1/4 .。。
int i = blockDim.x/2; while(i != 0) { if (cacheIndex < i) cache[cacheIndex] += cache[cacheIndex + i]; __syncthreads(); i /= 2; } if (cacheIndex == 0) c[blockIdx.x] = cache[0];
更好的优化算法：循环展开等
前缀和运算(Scan):
for(d = 1; (1 << d) < n; d++) for all k in parallel if( k >= (1 << d) ) x[out][k] = x[in][k – (1 << (d-1))] + x[in][k] else x[out][k] = x[in][k]
或
for(d = 1; (1 << d) < n; d++) for all k in parallel tmp = x[k] if( k >= (1 << d) ) tmp = x[k – (1 << (d-1))] + x[k] __syncthreads();//同步 x[k] = tmp
以上两算法运行所需空间至少是原空间的两倍，思想为倍增思想。
还有更高效的Scan算法。
for d:=0 to log2(n-1) do for k from 0 to n-1 by 2^(d+1) in parallel do x[k+2^(d+1)-1]:=x[k+2^(d+1)-1] + x[k+2^(d+1)-1] x[n-1]:=0 for d:=log2(n-1) downto 0 do for k from 0 to n-1 by 2^(d+1) in parallel do t:=x[k+2^d-1] x[k+2^d-1]:=x[k+2^(d+1)-1] x[k+2^(d+1)-1]:=t+x[k+2^(d+1)-1]
书上还有更高效的scan_best_kernel算法，略。
排序算法:
基于比较的排序：排序网络
基于非比较的排序：并行基数排序。前置技能：Scan。
并行基数排序算法： 按二进制位考虑。 以00101101为例。排完序后应当是12473568。 二进制翻转: 11010010 统计前缀和: 12233344 如果当前位是0，则写入对应位置。 第1个数写入首位置，第2个数写入第二个位置，第4个数写入第三个位置，第7个数写入第四个位置。 再对当前位是1的进行写入，位置下标 + 4(0的个数)。
矩阵乘法优化:
矩阵运算A*B = C， A为m*p的矩阵，B为p*n的矩阵。
优化1:
将C分块，每个线程块处理C中的某一块，例如d*d的一小块。
那么每个线程块需要完成d*p的矩阵与p*d的矩阵相乘的运算。
为了高效访存，每个线程块再对d*p和p*d的矩阵的p进行划分，看成多个矩阵块相乘后累加。
每个小块为d*q和q*d的大小，开在shared memory内，节约了大量global memory带宽。
（虽然循环次数会增加，但访存效率得到了高效提升）
优化2：
利用寄存器资源优化，效率更高，但略为繁琐。
矩阵转置优化:
无优化：拷贝至GPU内存，置换后拷贝回CPU内存。缺点：输入时每个block按行读入，满足合并访问条件；输出时数据间隔过大，不满足合并访问条件。
优化1：
分块，每个块是一个小方阵矩阵，如16*16。
输入时，每个线程块操作一个16*16方阵，通过shared memory完成16*16小方阵转置。
之后将大矩阵按块转置输出至global memory，每个线程块内无需再转置，满足合并访问条件。
shared memory数组大小设置成16*17而不是16*16，这样每行中处于同一列的数据就会被存储在不同的shared memory bank中，避免了bank conflict。
优化2：
上述无优化与优化1均存在分区冲突问题。优化2算法进行了for循环操作，暂未深入研究。
CUDA程序优化
grid和block的维度设计
grid的尺寸大一点较好。
为了有效利用执行单元，每个block的线程数应当是32的整数倍，最好让线程数量保持在64 ~ 256之间。
block维度和每个维度上的尺寸的主要作用是避免做整数除法和求模运算。实际中视具体情况而定。
如果问题规模对划分方式不敏感，应该让blockDim.x为16或16的整数倍，提高访问global memory和shared memory的效率。
存储器访问优化
灵活运用各种存储器的特性，实现最大可用带宽。
指令流优化
等
CUDA作业
作业1 简单CUDA应用，矩阵乘法
1 #include <bits/stdc++.h> 2 //#include "cuda_runtime.h" 3 //#include "device_launch_parameters.h" 4 5 using namespace std; 6 #define N 2000 7 8 const int block = 1<<12; 9 const int thread = 1<<10; 10 11 long long a[N][N]; 12 long long b[N][N]; 13 long long c[N][N]; 14 void init() { 15 for(int i = 0; i < N; i++) 16 for(int j = 0; j < N; j++) 17 a[i][j] = i*i+j, b[i][j] = i+j*j, c[i][j] = 0; 18 } 19 20 __global__ void init_cuda(long long *c) { 21 int id = blockIdx.x*blockDim.x+threadIdx.x; 22 if(id < N*N) c[id] = 0; 23 } 24 25 __global__ void mul_cuda(long long *a, long long *b, long long *c) { 26 int id = blockIdx.x*blockDim.x+threadIdx.x; 27 if(id < N*N) { 28 int row = id/N, col = id-row*N; 29 for(int k = 0; k < N; k++) 30 c[id] += a[row*N+k]*b[k*N+col]; 31 } 32 } 33 34 35 36 int main(int argc, char** argv) { 37 int cstart = clock(); 38 init(); 39 if(argv[1][0] == '0') { 40 puts("not cuda"); 41 for(int i = 0; i < N; i++) 42 for(int j = 0; j < N; j++) 43 for(int k = 0; k < N; k++) 44 c[i][k] += a[i][j]*b[j][k]; 45 } 46 else { 47 puts("cuda"); 48 long long *dev_a, *dev_b, *dev_c; 49 cudaMalloc( (void**)&dev_a, sizeof a ); 50 cudaMemcpy(dev_a, a, sizeof a, cudaMemcpyHostToDevice); 51 52 cudaMalloc( (void**)&dev_b, sizeof b ); 53 cudaMemcpy(dev_b, b, sizeof b, cudaMemcpyHostToDevice); 54 55 cudaMalloc( (void**)&dev_c, sizeof c ); 56 57 58 init_cuda<<<block, thread>>>(dev_c); 59 mul_cuda<<<block, thread>>>(dev_a, dev_b, dev_c); 60 61 cudaMemcpy(c, dev_c, sizeof c, cudaMemcpyDeviceToHost); 62 cudaFree(dev_a); 63 cudaFree(dev_b); 64 cudaFree(dev_c); 65 } 66 printf("%lld, ", c[1233][1233]); 67 printf("time: %d\n", int(clock()-cstart)); 68 return 0; 69 }
作业2 卷积操作，常量内存
1 //compile command: nvcc cv.cu `pkg-config --cflags --libs opencv` -std=c++11 2 //execute command1: ./a.out CC.jpg 3 3 //execute command2: ./a.out CC.jpg 5 4 #include <bits/stdc++.h> 5 6 #include <opencv2/opencv.hpp> 7 //#include <opencv2/gpu/gpu.hpp> 8 using namespace cv; 9 10 Mat G3 = (Mat_<int>(3, 3) << -5, 3, -5, 11 3, 9, 3, 12 -5, 3, -5); 13 14 Mat G5 = (Mat_<int>(5, 5) << 0, -1, 1, -1, 0, 15 -1, 1, -1, 1,-1, 16 0, -1, 8,-1, 0, 17 -1, 1, -1, 1,-1, 18 0, -1, 1,-1, 0); 19 20 void CPU_Sharpen(const Mat& myImage, Mat& Result, int ca){ 21 CV_Assert(myImage.depth() == CV_8U); // accept only uchar images 22 23 int begin = clock(); 24 Result.create(myImage.size(), myImage.type()); 25 const int nChannels = myImage.channels(); 26 Mat &G = ca == 3? G3: G5; 27 int half = G.rows >> 1; 28 29 for(int row = half; row < myImage.rows-half; ++row) { 30 uchar* output = Result.ptr<uchar>(row); 31 for(int col = half*nChannels; col < nChannels * (myImage.cols - half); ++col) { 32 int tmp = 0; 33 for(int i = 0; i < G.rows; ++i) 34 for(int j = 0; j < G.cols; ++j) 35 tmp += G.at<int>(i, j)*( *(myImage.ptr<uchar>(row-half+i)+(col-half*nChannels+j*nChannels) ) ); 36 *output++ = saturate_cast<uchar>(tmp); 37 } 38 } 39 for(int i = 0; i < half; i++) { 40 Result.row(i).setTo(Scalar(140)); 41 Result.row(Result.rows - 1 - i).setTo(Scalar(140)); 42 Result.col(i).setTo(Scalar(140)); 43 Result.col(Result.cols - 1 - i).setTo(Scalar(140)); 44 } 45 printf("Time used %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 46 } 47 48 /********************************************/ 49 50 __constant__ int con_G3[3][3] = { 51 {-5, 3, -5}, 52 { 3, 9, 3}, 53 {-5, 3, -5} 54 }; 55 __constant__ int con_G5[5][5] = { 56 {0, -1, 1,-1, 0}, 57 {-1, 1, -1, 1,-1}, 58 {0, -1, 8,-1, 0}, 59 {-1, 1, -1, 1,-1}, 60 {0, -1, 1,-1, 0} 61 }; 62 63 __global__ void init_cuda(uchar *c, int col_num) { 64 int col_id = blockIdx.x, row_id = threadIdx.x; 65 int now = (row_id*col_num+col_id)*3; 66 c[now] = c[now+1] = c[now+2] = 0; 67 } 68 69 70 //GPU,start from c, num * sizeof(uchar) 71 __global__ void test(uchar *c, int *sum, int num) { 72 int x = 0; 73 for(int i = 0; i < num; i++) 74 x += c[i]; 75 *sum = x; 76 } 77 78 __global__ void con_cuda(uchar *s, uchar *t, int ca, int row_num, int col_num) { 79 int col_id = blockIdx.x-1, row_id = threadIdx.x-1; 80 const int half = ca >> 1; 81 if(row_id >= half && row_id < row_num-half && col_id >= half && col_id < col_num-half) { 82 const int* con_mat = ca == 3? con_G3[0]: con_G5[0]; 83 int res[3] = {0, 0, 0}; 84 for(int i = 0; i < ca; i++) 85 for(int j = 0; j < ca; j++) { 86 //s[row_num][col_num][3]; 87 int pos = (row_id-half+i)*col_num*3+(col_id-half+j)*3; 88 res[0] += con_mat[i*ca+j]*s[pos]; 89 res[1] += con_mat[i*ca+j]*s[pos+1]; 90 res[2] += con_mat[i*ca+j]*s[pos+2]; 91 } 92 res[0] = res[0] < 0? 0: (res[0] > 255? 255: res[0]); 93 res[1] = res[1] < 0? 0: (res[1] > 255? 255: res[1]); 94 res[2] = res[2] < 0? 0: (res[2] > 255? 255: res[2]); 95 int pos = row_id*col_num*3+col_id*3; 96 t[pos] = res[0], 97 t[pos+1] = res[1], 98 t[pos+2] = res[2]; 99 } 100 } 101 102 /*******************************************/ 103 104 void HANDLE(cudaError x) { 105 if(x != cudaSuccess) { 106 puts("error!"); 107 exit(0); 108 } 109 } 110 111 int main(int argc, char** argv ) { 112 if ( argc < 3 ) { 113 printf("usage: a.out <Image_Path> <size of Mat>\n"); 114 return -1; 115 } 116 117 Mat src_img = imread(argv[1], 1), ans_CPU, ans_GPU; 118 int ca = argv[2][0]-'0'; 119 printf("%d %d\n", src_img.rows, src_img.cols); 120 121 /**********************************************************************************/ 122 123 printf("Run on CPU!\n"); 124 CPU_Sharpen(src_img, ans_CPU, ca); 125 std::string s = std::string("CC")+std::to_string(ca)+std::string("_With_CPU.jpg"); 126 imwrite(s, ans_CPU); 127 imshow("after operation", ans_CPU); 128 waitKey(); 129 130 /**********************************************************************************/ 131 132 printf("Run on GPU!\n"); 133 int begin = clock(); 134 uchar *dev_src, *dev_result; 135 int seg = src_img.cols*src_img.channels(); 136 HANDLE(cudaMalloc( (void**)&dev_src, src_img.rows*seg*sizeof(uchar))); 137 HANDLE(cudaMalloc( (void**)&dev_result, src_img.rows*seg*sizeof(uchar))); 138 /*Memcpy to dev_src*/ 139 for(int i = 0; i < src_img.rows; ++i) 140 HANDLE(cudaMemcpy(dev_src+i*seg*sizeof(uchar), src_img.ptr<uchar>(i), sizeof(uchar)*seg, cudaMemcpyHostToDevice)); 141 /*Init for dev_result*/ 142 init_cuda<<<src_img.cols, src_img.rows>>>(dev_result, src_img.cols); 143 /*Do convolution*/ 144 con_cuda<<<src_img.cols, src_img.rows>>>(dev_src, dev_result, ca, src_img.rows, src_img.cols); 145 146 ans_GPU.create(src_img.size(), src_img.type()); 147 /*Memcpy to host*/ 148 for(int i = 0; i < ans_GPU.rows; ++i) 149 cudaMemcpy(ans_GPU.ptr<uchar>(i), dev_result+i*seg*sizeof(uchar), sizeof(uchar)*seg, cudaMemcpyDeviceToHost); 150 151 for(int i = 0; i < (ca >> 1); i++) { 152 ans_GPU.row(i).setTo(Scalar(140)); 153 ans_GPU.row(ans_GPU.rows - 1 - i).setTo(Scalar(140)); 154 ans_GPU.col(i).setTo(Scalar(140)); 155 ans_GPU.col(ans_GPU.cols - 1 - i).setTo(Scalar(140)); 156 } 157 /*Free*/ 158 cudaFree(dev_src); 159 cudaFree(dev_result); 160 printf("Time used %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 161 imshow("after operation", ans_GPU); 162 waitKey(); 163 return 0; 164 }
作业3 卷积操作，流
1 //compile command: nvcc cv.cu `pkg-config --cflags --libs opencv` -std=c++11 2 //execute command1: ./a.out 1.jpg 3 3 //execute command2: ./a.out 1.jpg 5 4 #include <bits/stdc++.h> 5 6 #include <opencv2/opencv.hpp> 7 //#include <opencv2/gpu/gpu.hpp> 8 using namespace cv; 9 10 /********************************************/ 11 12 __constant__ int con_G3[3][3] = { 13 {-5, 3, -5}, 14 { 3, 9, 3}, 15 {-5, 3, -5} 16 }; 17 __constant__ int con_G5[5][5] = { 18 {0, -1, 1,-1, 0}, 19 {-1, 1, -1, 1,-1}, 20 {0, -1, 8,-1, 0}, 21 {-1, 1, -1, 1,-1}, 22 {0, -1, 1,-1, 0} 23 }; 24 25 __global__ void init_cuda(uchar *c, int col_num) { 26 int col_id = blockIdx.x, row_id = threadIdx.x; 27 int now = (row_id*col_num+col_id)*3; 28 c[now] = c[now+1] = c[now+2] = 0; 29 } 30 31 32 //GPU,start from c, num * sizeof(uchar) 33 __global__ void test(uchar *c, int *sum, int num) { 34 int x = 0; 35 for(int i = 0; i < num; i++) 36 x += c[i]; 37 *sum = x; 38 } 39 40 __global__ void con_cuda(uchar *s, uchar *t, int ca, int row_num, int col_num) { 41 int col_id = blockIdx.x-1, row_id = threadIdx.x-1; 42 const int half = ca >> 1; 43 if(row_id >= half && row_id < row_num-half && col_id >= half && col_id < col_num-half) { 44 const int* con_mat = ca == 3? con_G3[0]: con_G5[0]; 45 int res[3] = {0, 0, 0}; 46 for(int i = 0; i < ca; i++) 47 for(int j = 0; j < ca; j++) { 48 //s[row_num][col_num][3]; 49 int pos = (row_id-half+i)*col_num*3+(col_id-half+j)*3; 50 res[0] += con_mat[i*ca+j]*s[pos]; 51 res[1] += con_mat[i*ca+j]*s[pos+1]; 52 res[2] += con_mat[i*ca+j]*s[pos+2]; 53 } 54 res[0] = res[0] < 0? 0: (res[0] > 255? 255: res[0]); 55 res[1] = res[1] < 0? 0: (res[1] > 255? 255: res[1]); 56 res[2] = res[2] < 0? 0: (res[2] > 255? 255: res[2]); 57 int pos = row_id*col_num*3+col_id*3; 58 t[pos] = res[0], 59 t[pos+1] = res[1], 60 t[pos+2] = res[2]; 61 } 62 } 63 64 /*******************************************/ 65 66 void HANDLE_ERROR(cudaError x) { 67 if(x != cudaSuccess) { 68 puts("error!"); 69 exit(0); 70 } 71 } 72 73 int main(int argc, char** argv ) { 74 if ( argc < 3 ) { 75 printf("usage: a.out <Image_Path> <size of Mat>\n"); 76 return -1; 77 } 78 79 Mat src_img = imread(argv[1], 1), ans_CPU, ans_GPU; 80 int ca = argv[2][0]-'0'; 81 printf("%d %d\n", src_img.rows, src_img.cols); 82 /**********************************************************************************/ 83 84 printf("Run on GPU!\n"); 85 int begin = clock(); 86 /**********************************************************************************/ 87 88 uchar *dev_src, *dev_result; 89 int seg = src_img.cols*src_img.channels(); 90 HANDLE_ERROR(cudaMalloc( (void**)&dev_src, src_img.rows*seg*sizeof(uchar))); 91 HANDLE_ERROR(cudaMalloc( (void**)&dev_result, src_img.rows*seg*sizeof(uchar))); 92 93 cudaStream_t stream0, stream1; 94 HANDLE_ERROR( cudaStreamCreate( &stream0 ) ); 95 HANDLE_ERROR( cudaStreamCreate( &stream1 ) ); 96 for(int i = 0; i < src_img.rows; ++i) 97 HANDLE_ERROR( cudaMemcpyAsync(dev_src+i*seg*sizeof(uchar), src_img.ptr<uchar>(i), sizeof(uchar)*seg, cudaMemcpyHostToDevice, stream0) ); 98 99 init_cuda<<<src_img.cols, src_img.rows, 0, stream0>>>(dev_result, src_img.cols); 100 con_cuda<<<src_img.cols, src_img.rows, 0, stream0>>>(dev_src, dev_result, ca, src_img.rows, src_img.cols); 101 102 ans_GPU.create(src_img.size(), src_img.type()); 103 104 for(int i = 0; i < ans_GPU.rows; ++i) 105 HANDLE_ERROR( cudaMemcpyAsync(ans_GPU.ptr<uchar>(i), dev_result+i*seg*sizeof(uchar), sizeof(uchar)*seg, cudaMemcpyDeviceToHost, stream0) ); 106 107 for(int i = 0; i < (ca >> 1); i++) { 108 ans_GPU.row(i).setTo(Scalar(140)); 109 ans_GPU.row(ans_GPU.rows - 1 - i).setTo(Scalar(140)); 110 ans_GPU.col(i).setTo(Scalar(140)); 111 ans_GPU.col(ans_GPU.cols - 1 - i).setTo(Scalar(140)); 112 } 113 114 HANDLE_ERROR( cudaStreamSynchronize( stream0 ) ); 115 HANDLE_ERROR( cudaStreamSynchronize( stream1 ) ); 116 117 HANDLE_ERROR( cudaStreamDestroy( stream0 ) ); 118 HANDLE_ERROR( cudaStreamDestroy( stream1 ) ); 119 120 cudaFree(dev_src); 121 cudaFree(dev_result); 122 /********************without stream*********************/ 123 124 // uchar *dev_src, *dev_result; 125 // int seg = src_img.cols*src_img.channels(); 126 // HANDLE_ERROR(cudaMalloc( (void**)&dev_src, src_img.rows*seg*sizeof(uchar))); 127 // HANDLE_ERROR(cudaMalloc( (void**)&dev_result, src_img.rows*seg*sizeof(uchar))); 128 // /*Memcpy to dev_src*/ 129 // for(int i = 0; i < src_img.rows; ++i) 130 // HANDLE_ERROR(cudaMemcpy(dev_src+i*seg*sizeof(uchar), src_img.ptr<uchar>(i), sizeof(uchar)*seg, cudaMemcpyHostToDevice)); 131 // /*Init for dev_result*/ 132 // init_cuda<<<src_img.cols, src_img.rows>>>(dev_result, src_img.cols); 133 // /*Do convolution*/ 134 // con_cuda<<<src_img.cols, src_img.rows>>>(dev_src, dev_result, ca, src_img.rows, src_img.cols); 135 // 136 // ans_GPU.create(src_img.size(), src_img.type()); 137 // /*Memcpy to host*/ 138 // for(int i = 0; i < ans_GPU.rows; ++i) 139 // HANDLE_ERROR( cudaMemcpy(ans_GPU.ptr<uchar>(i), dev_result+i*seg*sizeof(uchar), sizeof(uchar)*seg, cudaMemcpyDeviceToHost) ); 140 // 141 // for(int i = 0; i < (ca >> 1); i++) { 142 // ans_GPU.row(i).setTo(Scalar(140)); 143 // ans_GPU.row(ans_GPU.rows - 1 - i).setTo(Scalar(140)); 144 // ans_GPU.col(i).setTo(Scalar(140)); 145 // ans_GPU.col(ans_GPU.cols - 1 - i).setTo(Scalar(140)); 146 // } 147 148 /*Free*/ 149 // cudaFree(dev_src); 150 // cudaFree(dev_result); 151 printf("Time used %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 152 imshow("after operation", ans_GPU); 153 imwrite("Tigerwith5.jpg", ans_GPU); 154 waitKey(); 155 return 0; 156 }
final project 图墙 + 图片融合
1 //compile command: nvcc final.cu `pkg-config --cflags --libs opencv` -std=c++11 2 //execute command1: ./a.out CC.jpg 3 3 //execute command2: ./a.out CC.jpg 5 4 #include <bits/stdc++.h> 5 6 #include <opencv2/opencv.hpp> 7 //#include <opencv2/gpu/gpu.hpp> 8 using namespace cv; 9 10 __global__ void init_cuda(uchar *c, int col_num) { 11 int col_id = blockIdx.x, row_id = threadIdx.x; 12 int now = (row_id*col_num+col_id)*3; 13 c[now] = c[now+1] = c[now+2] = 0; 14 } 15 16 17 //GPU,start from c, num * sizeof(uchar) 18 __global__ void test(uchar *c, int *sum, int num) { 19 int x = 0; 20 for(int i = 0; i < num; i++) 21 x += c[i]; 22 *sum = x; 23 } 24 25 __global__ void solve(uchar *s, uchar *t, uchar a, uchar b, uchar c, int R, int C) { 26 int x = blockIdx.x-1, y = threadIdx.x-1; 27 if(x < R&&y < C) { 28 int g = 3*(x*C+y); 29 if(t[g] == a&&t[g+1] == b&&t[g+2] == c) { 30 t[g] = s[g]; 31 t[g+1] = s[g+1]; 32 t[g+2] = s[g+2]; 33 } 34 else { 35 t[g] = 0.35*s[g] +0.65*t[g]; 36 t[g+1] = 0.35*s[g+1]+0.65*t[g+1]; 37 t[g+2] = 0.35*s[g+2]+0.65*t[g+2]; 38 } 39 } 40 } 41 42 __global__ void zoom(uchar *s, uchar *t, int R, int C, int r, int c) { 43 int x = blockIdx.x-1, y = threadIdx.x-1; 44 if(x <= r && y <= c) { 45 int row = x/(float)r*R, col = y/(float)c*C; 46 //t[x][y] = s[row][col]; 47 t[ (x*c+y)*3 ] = s[ (row*C+col)*3 ]; 48 t[ (x*c+y)*3+1 ] = s[ (row*C+col)*3+1 ]; 49 t[ (x*c+y)*3+2 ] = s[ (row*C+col)*3+2 ]; 50 } 51 } 52 53 /*******************************************/ 54 55 void HANDLE(cudaError x) { 56 if(x != cudaSuccess) { 57 puts("error!"); 58 exit(0); 59 } 60 } 61 62 const int N = 36; 63 64 int main(int argc, char** argv ) { 65 Mat src_img[N], dst_img[N], ret1, ret2; 66 int width = 1000000, height = 1000000; 67 for(int i = 0; i < N; i++) { 68 src_img[i] = imread(std::to_string(i+1)+std::string(".jpg"), 1); 69 int r = src_img[i].rows, c = src_img[i].cols; 70 if(height > r) height = r; 71 if(width > c) width = c; 72 } 73 height *= 0.4; 74 width *= 0.17; 75 76 //resize 77 int begin = clock(); 78 for(int i = 0; i < N; i++) { 79 dst_img[i].create(Size(height, width), src_img[i].type()); 80 resize(src_img[i], dst_img[i], Size(height, width)); 81 } 82 printf("Time used in resizing is%.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 83 84 int sq = sqrt(N+0.5); 85 //std::cout << width << ' ' << height << std::endl; 86 ret1.create(Size(height*sq, width*sq), src_img[0].type()); 87 //std::cout << ret1.rows << ' ' << ret1.cols << std::endl; 88 89 //merge 90 for(int i = 0; i < sq; i++) 91 for(int j = 0; j < sq; j++){ 92 for(int r = 0; r < width; r++) 93 memcpy(ret1.ptr<uchar>(i*width+r)+j*height*3, dst_img[i*sq+j].ptr<uchar>(r), height*3); 94 } 95 96 97 Mat ret = imread("0.jpg", 1); 98 resize(ret, ret2, Size(height*6, width*6)); 99 //std::cout << ret2.rows << ' ' << ret2.cols << std::endl; 100 //imshow("", ret2); 101 //waitKey(); 102 103 uchar a = *ret2.ptr<uchar>(0), b = *(ret2.ptr<uchar>(0)+1), c = *(ret2.ptr<uchar>(0)+2); 104 if(ret1.rows != ret2.rows || ret1.cols != ret2.cols) puts("gg"); 105 int R = ret2.rows, C = ret2.cols; 106 std::cout << R << ' ' << C << std::endl; 107 //CPU 108 begin = clock(); 109 for(int i = 0; i < R; i++) { 110 uchar *p1 = ret1.ptr<uchar>(i), *p2 = ret2.ptr<uchar>(i); 111 bool tag = true; 112 double x = 0; 113 for(int j = 0; j < C; j++) { 114 if(*(p2+j*3) == a&&*(p2+j*3+1) == b&&*(p2+j*3+2) == c) { 115 x = 0; 116 *(p2+j*3) = *(p1+j*3); 117 *(p2+j*3+1) = *(p1+j*3+1); 118 *(p2+j*3+2) = *(p1+j*3+2); 119 continue ; 120 } 121 122 if(*(p2+j*3+15) == a&&*(p2+j*3+16) == b&&*(p2+j*3+17) == c) { 123 x = 0; 124 *(p2+j*3) = *(p1+j*3); 125 *(p2+j*3+1) = *(p1+j*3+1); 126 *(p2+j*3+2) = *(p1+j*3+2); 127 continue ; 128 } 129 130 x = tag? x+0.06: x-0.06; 131 if(x > 1) tag = false; 132 if(x < 0.6) tag = true; 133 *(p2+j*3) = (1-x)*(*(p1+j*3)) +x*(*(p2+j*3)); 134 *(p2+j*3+1) = (1-x)*(*(p1+j*3+1))+x*(*(p2+j*3+1)); 135 *(p2+j*3+2) = (1-x)*(*(p1+j*3+2))+x*(*(p2+j*3+2)); 136 } 137 } 138 printf("Time used in CPU is %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 139 140 imshow("", ret2); 141 waitKey(); 142 imwrite("final.jpg", ret2); 143 //GPU 144 begin = clock(); 145 uchar *dev_src, *dev_result; 146 HANDLE(cudaMalloc( (void**)&dev_src, R*C*3*sizeof(uchar))); 147 HANDLE(cudaMalloc( (void**)&dev_result, R*C*3*sizeof(uchar))); 148 for(int i = 0; i < R; ++i) 149 cudaMemcpy(dev_src+i*C*3*sizeof(uchar), ret1.ptr<uchar>(i), C*3*sizeof(uchar), cudaMemcpyHostToDevice); 150 for(int i = 0; i < R; ++i) 151 cudaMemcpy(dev_result+i*C*3*sizeof(uchar), ret2.ptr<uchar>(i), C*3*sizeof(uchar), cudaMemcpyHostToDevice); 152 solve<<<R, C>>>(dev_src, dev_result, a, b, c, R, C); 153 for(int i = 0; i < R; ++i) 154 cudaMemcpy(ret2.ptr<uchar>(i), dev_result+i*C*3*sizeof(uchar), C*3*sizeof(uchar), cudaMemcpyDeviceToHost); 155 printf("Time used in GPU is %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 156 157 imshow("", ret2); 158 waitKey(); 159 cudaFree(dev_src); 160 cudaFree(dev_result); 161 imwrite("final2.jpg", ret2); 162 163 /**********************************************************************************/ 164 165 printf("Run on CPU!\n"); 166 CPU_Sharpen(src_img, ans_CPU, ca); 167 std::string s = std::string("IMG")+std::to_string(ca)+std::string("_With_CPU.jpg"); 168 imwrite(s, ans_CPU); 169 imshow("after operation", ans_CPU); 170 waitKey(); 171 172 /**********************************************************************************/ 173 174 printf("Run on GPU!\n"); 175 int begin = clock(); 176 uchar *dev_src, *dev_result; 177 int seg = src_img.cols*src_img.channels(); 178 HANDLE(cudaMalloc( (void**)&dev_src, src_img.rows*seg*sizeof(uchar))); 179 HANDLE(cudaMalloc( (void**)&dev_result, src_img.rows*seg*sizeof(uchar))); 180 181 /*Memcpy to dev_src*/ 182 183 for(int i = 0; i < src_img.rows; ++i) 184 HANDLE(cudaMemcpy(dev_src+i*seg*sizeof(uchar), src_img.ptr<uchar>(i), sizeof(uchar)*seg, cudaMemcpyHostToDevice)); 185 186 init_cuda<<<src_img.cols, src_img.rows>>>(dev_result, src_img.cols); 187 188 con_cuda<<<src_img.cols, src_img.rows>>>(dev_src, dev_result, ca, src_img.rows, src_img.cols); 189 190 ans_GPU.create(src_img.size(), src_img.type()); 191 192 for(int i = 0; i < ans_GPU.rows; ++i) 193 cudaMemcpy(ans_GPU.ptr<uchar>(i), dev_result+i*seg*sizeof(uchar), sizeof(uchar)*seg, cudaMemcpyDeviceToHost); 194 195 for(int i = 0; i < (ca >> 1); i++) { 196 ans_GPU.row(i).setTo(Scalar(140)); 197 ans_GPU.row(ans_GPU.rows - 1 - i).setTo(Scalar(140)); 198 ans_GPU.col(i).setTo(Scalar(140)); 199 ans_GPU.col(ans_GPU.cols - 1 - i).setTo(Scalar(140)); 200 } 201 202 cudaFree(dev_src); 203 cudaFree(dev_result); 204 printf("Time used %.3fms\n", ((int)clock()-begin)*1000.0/CLOCKS_PER_SEC); 205 imshow("after operation", ans_GPU); 206 waitKey(); 207 208 return 0; 209 }
emmmm，可能有bug，有冗余代码。
效果：
4万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
