DICOM形式の医用画像を用いて画像処理などを行いたい場合、当たり前ですがDICOM形式の医用画像ファイルを読まなければいけません。
世の中にはDICOM形式の医用画像の読み書きが出来るライブラリが諸々ありますが、例えばImageJでちょっとだけ加工してからコードで色々試したい、というようなときには、わざわざDICOM形式の画像をもう一度作ると面倒なので、値だけが直書きされたRawファイルでExportすることもそこそこあるかと思います(私だけ?)。
そして例えばCT画像の場合には各ピクセルのCT値が負の整数のときもありますので、Rawファイルは符号付16ビット整数で出力することになると思います。
というわけで、そんな符号付16ビット整数のバイナリファイルを読み書きするにはどうすれば良いのか、というお話です。
例として、4×4の以下のような値を各ピクセルに持つ画像を考えてみます。
画像内の数字は、各ピクセルの値で、普通の10進数表示と、対応する16進数表示です。
これをImageJでRawファイルとして出力してバイナリエディタで見てみると以下のようになっています。
ff c8 ff ed 01 01 07 06 00 e8 07 3f 0d 6a 10 cf 0c d3 11 20 10 91 0d a3 10 ed 0d dc 0c 4b 0b f6
16進数表示で順番に並んでいるだけの単純なものです。
これを効率よく読み込んで、配列に代入することを考えます。まずはPython版。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 |
import os import sys import array def read_signed_16bit_int_binary_sequence(file_full_path): file_size = os.stat(file_full_path).st_size # returns byte size number_of_pixels = file_size / 2 # 16 bit is 2 bytes. file = open(file_full_path, 'rb') image_value = array.array('h') # 'h': signed short (2 byte = 16 bit) image_value.fromfile(file, number_of_pixels) file.close() if sys.byteorder == 'little': image_value.byteswap() return image_value def write_signed_16bit_int_binary_sequence(file_full_path, value_array): file = open(file_full_path, 'wb') if sys.byteorder == 'little': value_array.byteswap() value_array.tofile(file) file.close() if sys.byteorder == 'little': value_array.byteswap() # read input_file_path = 'test.raw' pixel_values = read_signed_16bit_int_binary_sequence(input_file_path) # write output_file_path = 'test_output.raw' write_signed_16bit_int_binary_sequence(output_file_path, pixel_values) |
唯一の正解はありませんが、私が最も効率的だと思った方法はPython標準である
array.array('h')
を用いる方法です。これを用いると、
11 12 |
image_value = array.array('h') # 'h': signed short (2 byte = 16 bit) image_value.fromfile(file, number_of_pixels) |
の2行だけで、符号付16ビット整数が並んだバイナリファイルを一気に読み込むことが出来ます。
…だがしかし!その後の
15 16 |
if sys.byteorder == 'little': image_value.byteswap() |
この2行を忘れると、マシンによっては値がメチャクチャになります。Windowsマシンでは(※例外はあるかもしれません)具体的にはこうなってしまいます。
最初の図と16進数表示は同じなのに、対応する10進数表示がメチャクチャです。
array.array
を用いてバイナリファイルを読み込むと、そのバイナリファイルに書かれているそっくりそのままの状態がこのarray.arrayに書き込まれます。
いま、
11 |
image_value = array.array('h') # 'h': signed short (2 byte = 16 bit) |
と符号付short(Pythonではこれが符号付16ビット整数に対応)を指定していますが、これは2バイトなので、16進数2つずつが1つの整数として解釈されます。
image_value[0]
は、最初の16進数2つ、つまり
ff c8
から値が作られますが、メモリに配置されている通りに解釈されます。
ここでエンディアンの問題が発生します。Windowsはリトルエンディアンなので、メモリ上に
ff c8
と書かれている場合は、解釈としては
ffc8
ではなく、順番が逆になって
c8ff
になります。
なので、メモリ上に配置する際には各符号付16ビット整数で1バイト目と2バイト目とを入れ替えなければいけません。
プログラムを走らせているマシンがリトルエンディアンなのかビッグエンディアンなのかは、Pythonでは
sys.byteorder
という便利な機能で判定することが出来、バイトの入れ替えもarrayに対して
byteswap()
を行うと勝手にやってくれるという便利な機能があります。
次に書き込みですが、
array.array
の内容を全部バイナリファイルとして出力するための便利な機能がこれまた最初から用意されており、
27 |
value_array.tofile(file) |
の1行だけです。簡単ですね。
ただし、ここでもメモリの配置通りにファイルに書き込まれてしまいますので、もともとのファイルと同じにするためには、
byteswap()
を用いてエンディアンをビッグエンディアン表記にする必要があります。ファイルを書き込んだ後に再び
byteswap()
を行って、元に戻すのを忘れないようにしましょう。
Pythonではこのように、エンディアンの判定やbyteswapを簡単に行うことが出来て便利です。
が、Pythonはfor文がとても遅いため、画像全体に画像フィルタをかけたりするような場合には馬鹿にならない時間がかかります。
256×256の画像でも、ピクセル数は256×256 = 65,536もありますので、for文を回すと結構イライラしてきます。
Numbaなどを上手く用いれば高速化出来るのかもしれませんが…。
と言うわけで、C++でも書いてみました。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 |
#include<iostream> #include<fstream> #include<vector> #include<string> bool check_is_little_endian() { const int16_t value = 1; return (*(char*)&value) ? true : false; } const size_t get_file_byte_size(std::ifstream& file) { file.seekg(0, std::ios::end); const size_t file_size = (size_t)file.tellg(); file.seekg(0, std::ios::beg); return file_size; } void convert_signed_16bit_int_endian(std::vector<int16_t>* target_vector) { for (auto& value : *target_vector) { value = (value << 8) | ((value >> 8) & 0xFF); // The most left value does not change even if bit-shifted, so the mask by 0xFF is necessary. } } bool read_signed_16bit_int_binary_sequence(const std::string& file_full_path, std::vector<int16_t>* target_vector) { std::ifstream file(file_full_path, std::ios::in | std::ios::binary); if (!file) { std::cout << "Error: The file path is incorrect. There is no file." << std::endl; return false; } const size_t file_size = get_file_byte_size(file); target_vector->clear(); target_vector->resize(file_size / 2); // 16 bit is 2 bytes. file.read((char*)target_vector->data(), file_size); if (check_is_little_endian()) { convert_signed_16bit_int_endian(target_vector); } file.close(); return true; } bool write_signed_16bit_int_binary_sequence(const std::string& file_full_path, const std::vector<int16_t>& target_vector) { std::ofstream file(file_full_path, std::ios::out | std::ios::binary); if (!file) { std::cout << "Error: The file to write could not be opend." << std::endl; return false; } const size_t file_size = target_vector.size() * 2; // int16_t (16 bit) is 2 byte. if (check_is_little_endian()) { std::vector<int16_t> temp_vector = target_vector; // deep copy. convert_signed_16bit_int_endian(&temp_vector); file.write((char*)temp_vector.data(), file_size); } else { file.write((char*)target_vector.data(), file_size); } file.close(); return true; } int main() { const std::string input_file_path = "test.raw"; std::vector<int16_t> pixel_values; read_signed_16bit_int_binary_sequence(input_file_path, &pixel_values); const std::string output_file_path = "test_output.raw"; write_signed_16bit_int_binary_sequence(output_file_path, pixel_values); return 0; } |
Pythonと比べると行数が2倍弱も増えてしまっていますが、これは、Pythonでは標準で備わっていた機能のいくつかを自分で実装しなければいけないことも影響しています。
とは言えそこまで大変なことにはなりませんでした。
int16_t型の配列であるstd::vector<int16_t>に一気にバイナリファイルの内容を突っ込むためには
35 |
file.read((char*)target_vector->data(), file_size); |
の1行で事足りてしまいます。
ただ、ここでもやはりエンディアン問題が発生しますので、エンディアンを判定して、リトルエンディアンであればbyteswapが必要になります。
で、このエンディアン判定ですが、現在のC++の標準機能には無いようです。しかし判定が必要なときは少なからずありますので、誰かがうまいやり方を見つけているはずです。
検索すると簡単に見つかりました。
初心者のプログラミング体験記 ビッグ? or リトル?(エンディアンの判定)【C言語】
関数として実装すればこうなります。
6 7 8 9 |
bool check_is_little_endian() { const int16_t value = 1; return (*(char*)&value) ? true : false; } |
しかしbyteswapはちょっとハマりました…。
上位8ビットと下位8ビットとを入れ替えれば良いので、
21 |
value = (value << 8) | (value >> 8); |
でOK!…かと思いきや、これではダメでした。正しくは
21 |
value = (value << 8) | ((value >> 8) & 0xFF); // The most left value does not change even if bit-shifted, so the mask by 0xFF is necessary. |
です。Stack Overflowにありました。
& 0xFF
とは何ぞや?と思ったのですが、
整数型とビット操作
を読んで解決しました。符号付整数のビット演算では、ビットをずらしても符号の部分の値はそのまま残るので、上位8ビットを8ビット右にずらした場合は、その結果の16ビットのうちの一番左側のビットを常に0にしなければいけません。
そのための
& 0xFF
です。
そこさえクリアすれば、あとは書き込みの場合も含めてPythonと基本的には一緒です。
これでPythonでもC++でも符号付16ビット整数のバイナリファイルを自由に操作できますね。
もっと簡単な方法があればぜひ知りたいです!