DCM医学影像文件的格式与读取方法分享

转自 http://bbs.csdn.net/topics/390146622

2014-01-10 13:25 636人阅读 评论(0) 收藏 举报

医学影像DCM是个类似PNG的分块格式，内涵丰富医疗信息，由于应用领域较窄，OpenCV没有提供对其加载支持，一般采用dcmtk库进行加载。

博主本着研究探索的精神写了这段代码，并分享出来，提供了解析DCM的DIB数据、宽高、窗宽窗位、像素间距的功能，支持反色DCM。

将来可以参考spec文档进行扩充，可以在大框架中加入块的处理，从而抽取更多感兴趣的信息，也可以加入压缩DCM的支持。

[cpp] view plaincopy

1. // 输入：文件名  
2. // 输出：宽高 窗宽窗位 像素间距 dib  
3. unsigned short* dcmLoadImage(char* fn, int& width, int& height, int& windowWidth, int& windowLevel, double& pixelSpacing)  
4. {  
5.     // 读取整个文件到内存  
6.     FILE* f = fopen(fn, "rb");  
7.     fseek(f, 0, SEEK_END);  
8.     int dcmSize = ftell(f);  
9.     fseek(f, 0, SEEK_SET);  
10.     unsigned char* dcm = new unsigned char[dcmSize];  
11.     fread(dcm, dcmSize, 1, f);  
12.     fclose(f);  
13.       
14.     // 解析  
15.     width = -1;  
16.     height = -1;  
17.     windowWidth = -1;  
18.     windowLevel = -1;  
19.     pixelSpacing = 0.0;  
20.     bool invert = false;  
21.     unsigned short* dib = NULL;  
22.     int dibSize = -1;  
23.     int p = 132;  
24.     while( p<dcmSize )  
25.     {  
26.         short group = *(short*)(dcm+p);  
27.         p+=2;  
28.         short element = *(short*)(dcm+p);  
29.         p+=2;  
30.         //cout<<setfill('0')<<hex<<setw(4)<<group<<","<<setw(4)<<element<<setfill('\0')<<dec<<endl;  
31.   
32.   
33.         if(group==0x0002)  
34.         {  
35.             if(element==0x0001)// 2,1  
36.             {  
37.                 p+=10;  
38.             }  
39.             else// 2,*  
40.             {  
41.                 char vr[3];  
42.                 vr[0] = dcm[p++];  
43.                 vr[1] = dcm[p++];  
44.                 vr[2] = '\0';  
45.                 short size = *(short*)(dcm+p);  
46.                 p+=2;  
47.                 p+=size;  
48.             }  
49.         }  
50.         else// *  
51.         {  
52.             int size = *(int*)(dcm+p);  
53.             p+=4;  
54.             if(size==-1)  
55.                 size=0;  
56.             if(group==0x0028&&element==0x0010)  
57.                 height = *(short*)(dcm+p);  
58.             if(group==0x0028&&element==0x0030)  
59.                 pixelSpacing = atof((char*)dcm+p);  
60.             if(group==0x0028&&element==0x0004)  
61.                 invert = !strncmp((char*)dcm+p, "MONOCHROME1", 11);  
62.             else if(group==0x0028&&element==0x0011)  
63.                 width = *(short*)(dcm+p);  
64.             else if(group==0x0028&&element==0x1050)  
65.                 windowLevel = atoi((char*)dcm+p);  
66.             else if(group==0x0028&&element==0x1051)  
67.                 windowWidth = atoi((char*)dcm+p);  
68.             else if(group==0x7FE0&&element==0x0010)  
69.             {  
70.                 dibSize = size;  
71.                 dib = new unsigned short[dibSize];  
72.                 memcpy(dib, dcm+p, dibSize);  
73.             }  
74.             p+=size;  
75.         }  
76.     }  
77.     assert(p==dcmSize);  
78.     assert(width!=-1);  
79.     assert(height!=-1);  
80.     assert(windowLevel!=-1);  
81.     assert(windowWidth!=-1);  
82.     assert(dib!=NULL);  
83.     assert(dibSize==width*height*2);  
84.       
85.     // 反图做修正  
86.     if(invert)  
87.     {  
88.         for(int i=0; i<width*height; i++)  
89.             dib[i] ^= 0xFFFF;  
90.         windowLevel = 65535-windowLevel;  
91.     }  
92.       
93.     delete[] dcm;  
94.     return dib;  
95. }  
96.   
97.   
98. void dcmReleaseImage(unsigned short* dib)  
99. {  
100.     delete[] dib;  
101. }