1 /*=========================================================================
4 Module: $RCSfile: Dense2007ToDicom.cxx,v $
6 Date: $Date: 2008/03/17 13:16:10 $
7 Version: $Revision: 1.6 $
9 Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
10 l'Image). All rights reserved. See Doc/License.txt or
11 http://www.creatis.insa-lyon.fr/Public/Gdcm/License.html for details.
13 This software is distributed WITHOUT ANY WARRANTY; without even
14 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15 PURPOSE. See the above copyright notices for more information.
17 =========================================================================*/
23 #if defined(__BORLANDC__)
28 #include "gdcmFileHelper.h"
29 #include "gdcmDebug.h"
30 #include "gdcmDirList.h"
32 #include "gdcmArgMgr.h"
36 * Converts the "Dense" ".txt" (2007 version) files into 16 bits Dicom Files,
37 * Hope they don't change soon!
41 void LoadPeakStrain(std::ifstream &from, std::string imageName);
42 void LoadStrain(std::ifstream &from, std::string imageName);
43 void MakeDicomImage(float *tabVal, float *X, float *Y, float *Z, int NP, std::string dcmImageName);
47 int main(int argc, char *argv[])
50 " \n Dense2007ToDicom :\n ",
51 " Converts the '.txt' files into 16 bits Dicom Files, ",
53 " Dense2007ToDicom strain=...strain.txt peak_strain=...peak_strain.txt ",
54 " [verbose] [debug] ",
56 " verbose : user wants to run the program in 'verbose mode' ",
57 " debug : *developer* wants to run the program in 'debug mode' ",
60 // ----- Initialize Arguments Manager ------
62 GDCM_NAME_SPACE::ArgMgr *am = new GDCM_NAME_SPACE::ArgMgr(argc, argv);
64 if (argc == 1 || am->ArgMgrDefined("usage"))
66 am->ArgMgrUsage(usage); // Display 'usage'
70 // Seems that ArgMgrWantString doesn't work on MacOS
71 if(!am->ArgMgrDefined("strain"))
73 std::cout << "strain is mandatory" << std::endl;
76 if(!am->ArgMgrDefined("peak_strain"))
78 std::cout << "peak_strain is mandatory" << std::endl;
82 const char *strain = am->ArgMgrWantString("strain",usage);
83 const char *peak_strain = am->ArgMgrWantString("peak_strain",usage);
85 if (am->ArgMgrDefined("debug"))
86 GDCM_NAME_SPACE::Debug::DebugOn();
88 verbose = am->ArgMgrDefined("verbose");
90 // if unused Param we give up
91 if ( am->ArgMgrPrintUnusedLabels() )
93 am->ArgMgrUsage(usage);
97 delete am; // we don't need Argument Manager any longer
99 // ----- Begin Processing -----
102 std::ifstream fromPeakStrain( peak_strain );
103 if ( !fromPeakStrain )
105 std::cout << "Can't open file" << peak_strain << std::endl;
109 std::ifstream fromStrain( strain );
112 std::cout << "Can't open file" << strain << std::endl;
116 std::cout << "Success in open file" << peak_strain << std::endl;
117 LoadPeakStrain(fromPeakStrain, peak_strain);
118 fromPeakStrain.close();
120 std::cout << "Success in open file" << strain << std::endl;
121 LoadStrain(fromStrain, strain);
126 // =====================================================================================================================
128 void LoadPeakStrain(std::ifstream &from, std::string textFileName)
130 // in sax_base_slice0_peak_strain.txt :
133 Number of material points (NP) = 181
134 Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.3243 3.19392
136 Readout direction in 3D = -0.162314 -0.0771294 -0.983720
137 Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
138 Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
139 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
140 followed by their peak Ecc strain, an array of NP elements,
141 followed by their peak Err strain, an array of NP elements,
142 followed by their peak E11 strain, an array of NP elements,
143 followed by their Peak E22 strain, an array of NP elements,
144 42.0000 10.0000 0.000000
146 -0.154905 -0.0840482 -0.157350 -0.221403 -0.168118 -0.131331
147 -0.153781 -0.148481 -0.166602 -0.232858 -0.222650 -0.213712
157 //Number of material points (NP) = 181
166 std::cout << "NP : " << NP << std::endl;
168 //Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.3243 3.19392 88.2381
181 float readout, phase_enc, slice_sel;
185 std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
187 // Readout direction in 3D = -0.162314 -0.0771294 -0.983720
195 float readoutX, readoutY, readoutZ;
199 std::cout << " readoutX " << readoutX << " readoutY " << readoutY << " readoutZ " << readoutZ << std::endl;
201 // Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
210 float phase_encX, phase_encY, phase_encZ;
214 std::cout << " phase_encX " << phase_encX << " phase_encY " << phase_encY << " phase_encZ " << phase_encZ << std::endl;
216 // Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
224 float slice_selX, slice_selY, slice_selZ;
228 std::cout << " slice_selX " << slice_selX << " slice_selY " << slice_selY << " slice_selZ " << slice_selZ << std::endl;
234 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
235 followed by their peak Ecc strain, an array of NP elements,
236 followed by their peak Err strain, an array of NP elements,
237 followed by their peak E11 strain, an array of NP elements,
238 followed by their Peak E22 strain, an array of NP elements,
241 std::cout << "------------start skipping 1 line---------------- " << std::endl;
242 std::getline(from, str1);
243 std::cout << "[" << str1 << "]" << std::endl;
244 std::cout << "------------start skipping 1 line---------------- " << std::endl;
245 std::getline(from, str1);
246 std::cout << "[" << str1 << "]" << std::endl;
247 std::cout << "------------start skipping 1 line---------------- " << std::endl;
248 std::getline(from, str1);
249 std::cout << "[" << str1 << "]" << std::endl;
250 std::cout << "------------start skipping 1 line---------------- " << std::endl;
251 std::getline(from, str1);
252 std::cout << "[" << str1 << "]" << std::endl;
253 std::cout << "------------start skipping 1 line---------------- " << std::endl;
254 std::getline(from, str1);
255 std::cout << "[" << str1 << "]" << std::endl;
256 std::cout << "------------start skipping 1 line---------------- " << std::endl;
257 std::getline(from, str1);
258 std::cout << "[" << str1 << "]" << std::endl;
259 std::cout << "------------stop skipping ---------------- " << std::endl;
261 float *X = new float[NP];
262 float *Y = new float[NP];
263 float *Z = new float[NP];
267 for (i=0; i<NP; i++) {
292 std::cout << "--------------- Ecc_strain ------------------" << std::endl;
293 float *ecc_strain = new float[NP];
294 for (i=0; i<NP; i++) {
295 from >> ecc_strain[i];
297 std::cout << ecc_strain[i] << std::endl;
300 std::cout << "--------------- Err_strain ------------------" << std::endl;
301 float *err_strain = new float[NP];
302 for (i=0; i<NP; i++) {
303 from >> err_strain[i];
305 std::cout << err_strain[i] << std::endl;
308 std::cout << "--------------- E11_strain ------------------" << std::endl;
309 float *e11_strain = new float[NP];
310 for (i=0; i<NP; i++) {
311 from >> e11_strain[i];
313 std::cout << e11_strain[i] << std::endl;
316 std::cout << "--------------- E22_strain ------------------" << std::endl;
317 float *e22_strain = new float[NP];
318 for (i=0; i<NP; i++) {
319 from >> e22_strain[i];
321 std::cout << e22_strain[i] << std::endl;
324 std::string dcmImageName;
326 //followed by their peak Ecc strain, an array of NP elements,
327 dcmImageName = textFileName + "_peak_Ecc_strain.dcm";
328 MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
330 //followed by their peak Err strain, an array of NP elements,
331 dcmImageName = textFileName + "_peak_Err_strain.dcm";
332 MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
334 //followed by their peak E11 strain, an array of NP elements,
335 dcmImageName = textFileName + "_peak_E11_strain.dcm";
336 MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
338 //followed by their Peak E22 strain, an array of NP elements,
339 dcmImageName = textFileName + "_peak_E22_strain.dcm";
340 MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
343 // =====================================================================================================================
345 void LoadStrain(std::ifstream &from, std::string textFileName)
348 // in sax_base_slice0_strain.txt :
350 Number of cine frames = 18
351 Temporal resolution = 32.0000 ms
352 First frame starts at 48.0000 ms
353 Number of material points (NP) = 181
354 Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.324341 3.193918 88.238113
355 Readout direction in 3D = -0.162314 -0.0771294 -0.983720
356 Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
357 Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
358 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
359 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
360 followed by their Err strain, an array of dimensions(NP, number of cine frames),
361 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
362 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
363 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
364 42.0000 10.0000 0.000000
365 44.0000 10.0000 0.000000
367 -0.0622793 -0.0840482 -0.157350 -0.196722 -0.105844 -0.131331
368 -0.153781 -0.00940573 -0.0542236 -0.100403 -0.0369671 -0.0696840
375 int NP; // Number of Pints
376 int NCF; // Number of cine frames
377 float TR; // Temporal resolution
378 float FFS; // First frame starts
380 // Number of cine frames = 18
388 // Temporal resolution = 32.0000 ms
395 // First frame starts at 48.0000 ms
403 //Number of material points (NP) = 181
412 std::cout << "NP : " << NP << std::endl;
414 //Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.324341 3.193918 88.238113
427 float readout, phase_enc, slice_sel;
431 std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
433 // Readout direction in 3D = -0.162314 -0.0771294 -0.983720
441 float readoutX, readoutY, readoutZ;
445 std::cout << " readoutX " << readoutX << " readoutY " << readoutY << " readoutZ " << readoutZ << std::endl;
447 // Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
456 float phase_encX, phase_encY, phase_encZ;
460 std::cout << " phase_encX " << phase_encX << " phase_encY " << phase_encY << " phase_encZ " << phase_encZ << std::endl;
462 // Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
470 float slice_selX, slice_selY, slice_selZ;
474 std::cout << " slice_selX " << slice_selX << " slice_selY " << slice_selY << " slice_selZ " << slice_selZ << std::endl;
480 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
481 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
482 followed by their Err strain, an array of dimensions(NP, number of cine frames),
483 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
484 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
485 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
487 std::cout << "------------start skipping 1 line---------------- " << std::endl;
488 std::getline(from, str1);
489 std::cout << "[" << str1 << "]" << std::endl;
490 std::cout << "------------start skipping 1 line---------------- " << std::endl;
491 std::getline(from, str1);
492 std::cout << "[" << str1 << "]" << std::endl;
493 std::cout << "------------start skipping 1 line---------------- " << std::endl;
494 std::getline(from, str1);
495 std::cout << "[" << str1 << "]" << std::endl;
496 std::cout << "------------start skipping 1 line---------------- " << std::endl;
497 std::getline(from, str1);
498 std::cout << "[" << str1 << "]" << std::endl;
499 std::cout << "------------start skipping 1 line---------------- " << std::endl;
500 std::getline(from, str1);
501 std::cout << "[" << str1 << "]" << std::endl;
502 std::cout << "------------start skipping 1 line---------------- " << std::endl;
503 std::getline(from, str1);
504 std::cout << "[" << str1 << "]" << std::endl;
505 std::cout << "------------stop skipping ---------------- " << std::endl;
506 std::getline(from, str1);
507 std::cout << "[" << str1 << "]" << std::endl;
508 std::cout << "------------stop skipping ---------------- " << std::endl;
510 float *X = new float[NP];
511 float *Y = new float[NP];
512 float *Z = new float[NP];
516 for (i=0; i<NP; i++) {
539 std::cout << "--------------- Ecc_strain ------------------" << std::endl;
540 float *ecc_strain = new float[NP];
541 for (i=0; i<NP; i++) {
542 from >> ecc_strain[i];
544 std::cout << ecc_strain[i] << std::endl;
547 std::cout << "--------------- Err_strain ------------------" << std::endl;
548 float *err_strain = new float[NP];
549 for (i=0; i<NP; i++) {
550 from >> err_strain[i];
552 std::cout << err_strain[i] << std::endl;
555 std::cout << "--------------- E11_strain ------------------" << std::endl;
556 float *e11_strain = new float[NP];
557 for (i=0; i<NP; i++) {
558 from >> e11_strain[i];
560 std::cout << e11_strain[i] << std::endl;
563 std::cout << "--------------- E22_strain ------------------" << std::endl;
564 float *e22_strain = new float[NP];
565 for (i=0; i<NP; i++) {
566 from >> e22_strain[i];
568 std::cout << e22_strain[i] << std::endl;
571 std::string dcmImageName;
573 //followed by their Ecc strain, an array of NP elements,
574 dcmImageName = textFileName + "_Ecc_strain.dcm";
575 MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
577 //followed by their Err strain, an array of NP elements,
578 dcmImageName = textFileName + "_Err_strain.dcm";
579 MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
581 //followed by their E11 strain, an array of NP elements,
582 dcmImageName = textFileName + "_E11_strain.dcm";
583 MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
585 //followed by their E22 strain, an array of NP elements,
586 dcmImageName = textFileName + "_E22_strain.dcm";
587 MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
592 // =====================================================================================================================
595 void MakeDicomImage(float *tabVal, float *X, float *Y, float *Z, int NP, std::string dcmImageName)
597 float minX = 99999., minY = 99999., minZ = 99999.;
598 float maxX = 0., maxY = 0., maxZ = 0.;
601 for (i=0; i<NP; i++) {
602 // std::cout << X[i] << " " << Y[i] << " " << Z[i] << std::endl;
617 std::cout << "Min X,Y,Z " << minX << " " << minY << " " << minZ << std::endl;
618 std::cout << "Max X,Y,Z " << maxX << " " << maxY << " " << maxZ << std::endl;
619 std::cout << "Size X,Y,Z " << maxX-minX << " " << maxY-minY << " " << maxZ-minZ << std::endl;
621 // uint16_t *img = new uint16_t[int(maxX+0.5)*int(maxY+0.5)];
622 uint16_t *img = new uint16_t[int(maxX*4.)*int(maxY*4.)];
624 // Set whole image to 0
625 for(int i3=0;i3<int(maxX*4.)*int(maxY*4.);i3++)
626 // for(int i3=0;i3<int(maxX+0.5)*int(maxY+0.5);i3++)
629 for(int i2=0; i2<NP; i2++) {
631 int ordX = int(X[i2]*4.-30);
632 int ordY = int(maxY*4.) - int(Y[i2]*4.)+30;
633 // img[ /*int(maxX) -*/ int(X[i2]+0.5-1) + (int(maxY+0.5) - int(Y[i2]+0.5)) * int(maxX+0.5) ] = int(tabVal[i2]*100);
635 // img[ /*int(maxX) -*/ int(X[i2]*4.-30) + (int(maxY*4.) - int(Y[i2]*4.)+30) * int(maxX*4.) ] = int(tabVal[i2]*100);
636 img[ /*int(maxX) -*/ ordX + ordY * int(maxX*4.) ] = int(tabVal[i2]*100);
638 // Try to round up, just to see.
639 for(int iii=ordY-3; iii<ordY+4; iii++) for(int jjj=ordX-3; jjj<ordX+4; jjj++)
640 img[ jjj + iii * int(maxX*4.) ] = int(tabVal[i2]*100);
641 std::cout << int(X[i2]*4.) << " " << int(Y[i2]*4.) << " = " << int(tabVal[i2]*100) << std::endl;
644 // GDCM_NAME_SPACE::Debug::DebugOn();
646 std::ostringstream str;
648 GDCM_NAME_SPACE::File *file;
649 file = GDCM_NAME_SPACE::File::New();
651 // Set the image size
653 str << (int)(maxX*4.);
654 file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
656 str << (int)(maxY*4.);
657 file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
659 // Set the pixel type
661 file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
663 str << 16; // may be 12 or 16 if componentSize =16
664 file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
665 file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
667 // Set the pixel representation // 0/1 , 0=unsigned
668 file->InsertEntryString("1",0x0028,0x0103, "US"); // Pixel Representation
669 // Set the samples per pixel // 1:Grey level, 3:RGB
670 file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
673 // Set Rescale Intercept
676 file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
681 file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
684 GDCM_NAME_SPACE::FileHelper *fileH;
685 fileH = GDCM_NAME_SPACE::FileHelper::New(file);
686 // cast is just to avoid warnings (*no* conversion)
687 //fileH->SetImageData((uint8_t *)img,int(maxX*maxY)*sizeof(uint16_t)); // troubles when maxX, mayY are *actually* float!
688 fileH->SetImageData((uint8_t *)img,int(maxX*4.)*int(maxY*4.)*sizeof(uint16_t));
689 fileH->SetWriteModeToRaw();
690 fileH->SetWriteTypeToDcmExplVR();
692 if( !fileH->Write(dcmImageName))
693 std::cout << "Failed for [" << dcmImageName << "]\n"
694 << " File is unwrittable" << std::endl;