1 /*=========================================================================
4 Module: $RCSfile: Dense2007ToDicom.cxx,v $
6 Date: $Date: 2007/10/29 17:13:59 $
7 Version: $Revision: 1.5 $
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'
71 const char *strain = am->ArgMgrWantString("strain",usage);
72 const char *peak_strain = am->ArgMgrWantString("peak_strain",usage);
74 if (am->ArgMgrDefined("debug"))
75 GDCM_NAME_SPACE::Debug::DebugOn();
77 verbose = am->ArgMgrDefined("verbose");
79 // if unused Param we give up
80 if ( am->ArgMgrPrintUnusedLabels() )
82 am->ArgMgrUsage(usage);
86 delete am; // we don't need Argument Manager any longer
88 // ----- Begin Processing -----
91 std::ifstream fromPeakStrain( peak_strain );
92 if ( !fromPeakStrain )
94 std::cout << "Can't open file" << peak_strain << std::endl;
98 std::ifstream fromStrain( strain );
101 std::cout << "Can't open file" << strain << std::endl;
105 std::cout << "Success in open file" << peak_strain << std::endl;
106 LoadPeakStrain(fromPeakStrain, peak_strain);
107 fromPeakStrain.close();
109 std::cout << "Success in open file" << strain << std::endl;
110 LoadStrain(fromStrain, strain);
115 // =====================================================================================================================
117 void LoadPeakStrain(std::ifstream &from, std::string textFileName)
119 // in sax_base_slice0_peak_strain.txt :
122 Number of material points (NP) = 181
123 Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.3243 3.19392
125 Readout direction in 3D = -0.162314 -0.0771294 -0.983720
126 Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
127 Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
128 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
129 followed by their peak Ecc strain, an array of NP elements,
130 followed by their peak Err strain, an array of NP elements,
131 followed by their peak E11 strain, an array of NP elements,
132 followed by their Peak E22 strain, an array of NP elements,
133 42.0000 10.0000 0.000000
135 -0.154905 -0.0840482 -0.157350 -0.221403 -0.168118 -0.131331
136 -0.153781 -0.148481 -0.166602 -0.232858 -0.222650 -0.213712
146 //Number of material points (NP) = 181
155 std::cout << "NP : " << NP << std::endl;
157 //Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.3243 3.19392 88.2381
170 float readout, phase_enc, slice_sel;
174 std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
176 // Readout direction in 3D = -0.162314 -0.0771294 -0.983720
184 float readoutX, readoutY, readoutZ;
188 std::cout << " readoutX " << readoutX << " readoutY " << readoutY << " readoutZ " << readoutZ << std::endl;
190 // Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
199 float phase_encX, phase_encY, phase_encZ;
203 std::cout << " phase_encX " << phase_encX << " phase_encY " << phase_encY << " phase_encZ " << phase_encZ << std::endl;
205 // Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
213 float slice_selX, slice_selY, slice_selZ;
217 std::cout << " slice_selX " << slice_selX << " slice_selY " << slice_selY << " slice_selZ " << slice_selZ << std::endl;
223 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
224 followed by their peak Ecc strain, an array of NP elements,
225 followed by their peak Err strain, an array of NP elements,
226 followed by their peak E11 strain, an array of NP elements,
227 followed by their Peak E22 strain, an array of NP elements,
230 std::cout << "------------start skipping 1 line---------------- " << std::endl;
231 std::getline(from, str1);
232 std::cout << "[" << str1 << "]" << std::endl;
233 std::cout << "------------start skipping 1 line---------------- " << std::endl;
234 std::getline(from, str1);
235 std::cout << "[" << str1 << "]" << std::endl;
236 std::cout << "------------start skipping 1 line---------------- " << std::endl;
237 std::getline(from, str1);
238 std::cout << "[" << str1 << "]" << std::endl;
239 std::cout << "------------start skipping 1 line---------------- " << std::endl;
240 std::getline(from, str1);
241 std::cout << "[" << str1 << "]" << std::endl;
242 std::cout << "------------start skipping 1 line---------------- " << std::endl;
243 std::getline(from, str1);
244 std::cout << "[" << str1 << "]" << std::endl;
245 std::cout << "------------start skipping 1 line---------------- " << std::endl;
246 std::getline(from, str1);
247 std::cout << "[" << str1 << "]" << std::endl;
248 std::cout << "------------stop skipping ---------------- " << std::endl;
250 float *X = new float[NP];
251 float *Y = new float[NP];
252 float *Z = new float[NP];
256 for (i=0; i<NP; i++) {
281 std::cout << "--------------- Ecc_strain ------------------" << std::endl;
282 float *ecc_strain = new float[NP];
283 for (i=0; i<NP; i++) {
284 from >> ecc_strain[i];
286 std::cout << ecc_strain[i] << std::endl;
289 std::cout << "--------------- Err_strain ------------------" << std::endl;
290 float *err_strain = new float[NP];
291 for (i=0; i<NP; i++) {
292 from >> err_strain[i];
294 std::cout << err_strain[i] << std::endl;
297 std::cout << "--------------- E11_strain ------------------" << std::endl;
298 float *e11_strain = new float[NP];
299 for (i=0; i<NP; i++) {
300 from >> e11_strain[i];
302 std::cout << e11_strain[i] << std::endl;
305 std::cout << "--------------- E22_strain ------------------" << std::endl;
306 float *e22_strain = new float[NP];
307 for (i=0; i<NP; i++) {
308 from >> e22_strain[i];
310 std::cout << e22_strain[i] << std::endl;
313 std::string dcmImageName;
315 //followed by their peak Ecc strain, an array of NP elements,
316 dcmImageName = textFileName + "_peak_Ecc_strain.dcm";
317 MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
319 //followed by their peak Err strain, an array of NP elements,
320 dcmImageName = textFileName + "_peak_Err_strain.dcm";
321 MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
323 //followed by their peak E11 strain, an array of NP elements,
324 dcmImageName = textFileName + "_peak_E11_strain.dcm";
325 MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
327 //followed by their Peak E22 strain, an array of NP elements,
328 dcmImageName = textFileName + "_peak_E22_strain.dcm";
329 MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
332 // =====================================================================================================================
334 void LoadStrain(std::ifstream &from, std::string textFileName)
337 // in sax_base_slice0_strain.txt :
339 Number of cine frames = 18
340 Temporal resolution = 32.0000 ms
341 First frame starts at 48.0000 ms
342 Number of material points (NP) = 181
343 Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.324341 3.193918 88.238113
344 Readout direction in 3D = -0.162314 -0.0771294 -0.983720
345 Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
346 Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
347 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
348 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
349 followed by their Err strain, an array of dimensions(NP, number of cine frames),
350 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
351 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
352 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
353 42.0000 10.0000 0.000000
354 44.0000 10.0000 0.000000
356 -0.0622793 -0.0840482 -0.157350 -0.196722 -0.105844 -0.131331
357 -0.153781 -0.00940573 -0.0542236 -0.100403 -0.0369671 -0.0696840
364 int NP; // Number of Pints
365 int NCF; // Number of cine frames
366 float TR; // Temporal resolution
367 float FFS; // First frame starts
369 // Number of cine frames = 18
377 // Temporal resolution = 32.0000 ms
384 // First frame starts at 48.0000 ms
392 //Number of material points (NP) = 181
401 std::cout << "NP : " << NP << std::endl;
403 //Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.324341 3.193918 88.238113
416 float readout, phase_enc, slice_sel;
420 std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
422 // Readout direction in 3D = -0.162314 -0.0771294 -0.983720
430 float readoutX, readoutY, readoutZ;
434 std::cout << " readoutX " << readoutX << " readoutY " << readoutY << " readoutZ " << readoutZ << std::endl;
436 // Phase Enc. direction in 3D = -0.540606 -0.827052 0.154046
445 float phase_encX, phase_encY, phase_encZ;
449 std::cout << " phase_encX " << phase_encX << " phase_encY " << phase_encY << " phase_encZ " << phase_encZ << std::endl;
451 // Slice select direction in 3D = 0.825469 -0.556809 -0.0925458
459 float slice_selX, slice_selY, slice_selZ;
463 std::cout << " slice_selX " << slice_selX << " slice_selY " << slice_selY << " slice_selZ " << slice_selZ << std::endl;
469 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
470 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
471 followed by their Err strain, an array of dimensions(NP, number of cine frames),
472 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
473 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
474 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
476 std::cout << "------------start skipping 1 line---------------- " << std::endl;
477 std::getline(from, str1);
478 std::cout << "[" << str1 << "]" << std::endl;
479 std::cout << "------------start skipping 1 line---------------- " << std::endl;
480 std::getline(from, str1);
481 std::cout << "[" << str1 << "]" << std::endl;
482 std::cout << "------------start skipping 1 line---------------- " << std::endl;
483 std::getline(from, str1);
484 std::cout << "[" << str1 << "]" << std::endl;
485 std::cout << "------------start skipping 1 line---------------- " << std::endl;
486 std::getline(from, str1);
487 std::cout << "[" << str1 << "]" << std::endl;
488 std::cout << "------------start skipping 1 line---------------- " << std::endl;
489 std::getline(from, str1);
490 std::cout << "[" << str1 << "]" << std::endl;
491 std::cout << "------------start skipping 1 line---------------- " << std::endl;
492 std::getline(from, str1);
493 std::cout << "[" << str1 << "]" << std::endl;
494 std::cout << "------------stop skipping ---------------- " << std::endl;
495 std::getline(from, str1);
496 std::cout << "[" << str1 << "]" << std::endl;
497 std::cout << "------------stop skipping ---------------- " << std::endl;
499 float *X = new float[NP];
500 float *Y = new float[NP];
501 float *Z = new float[NP];
505 for (i=0; i<NP; i++) {
528 std::cout << "--------------- Ecc_strain ------------------" << std::endl;
529 float *ecc_strain = new float[NP];
530 for (i=0; i<NP; i++) {
531 from >> ecc_strain[i];
533 std::cout << ecc_strain[i] << std::endl;
536 std::cout << "--------------- Err_strain ------------------" << std::endl;
537 float *err_strain = new float[NP];
538 for (i=0; i<NP; i++) {
539 from >> err_strain[i];
541 std::cout << err_strain[i] << std::endl;
544 std::cout << "--------------- E11_strain ------------------" << std::endl;
545 float *e11_strain = new float[NP];
546 for (i=0; i<NP; i++) {
547 from >> e11_strain[i];
549 std::cout << e11_strain[i] << std::endl;
552 std::cout << "--------------- E22_strain ------------------" << std::endl;
553 float *e22_strain = new float[NP];
554 for (i=0; i<NP; i++) {
555 from >> e22_strain[i];
557 std::cout << e22_strain[i] << std::endl;
560 std::string dcmImageName;
562 //followed by their Ecc strain, an array of NP elements,
563 dcmImageName = textFileName + "_Ecc_strain.dcm";
564 MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
566 //followed by their Err strain, an array of NP elements,
567 dcmImageName = textFileName + "_Err_strain.dcm";
568 MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
570 //followed by their E11 strain, an array of NP elements,
571 dcmImageName = textFileName + "_E11_strain.dcm";
572 MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
574 //followed by their E22 strain, an array of NP elements,
575 dcmImageName = textFileName + "_E22_strain.dcm";
576 MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
581 // =====================================================================================================================
584 void MakeDicomImage(float *tabVal, float *X, float *Y, float *Z, int NP, std::string dcmImageName)
586 float minX = 99999., minY = 99999., minZ = 99999.;
587 float maxX = 0., maxY = 0., maxZ = 0.;
590 for (i=0; i<NP; i++) {
591 // std::cout << X[i] << " " << Y[i] << " " << Z[i] << std::endl;
606 std::cout << "Min X,Y,Z " << minX << " " << minY << " " << minZ << std::endl;
607 std::cout << "Max X,Y,Z " << maxX << " " << maxY << " " << maxZ << std::endl;
608 std::cout << "Size X,Y,Z " << maxX-minX << " " << maxY-minY << " " << maxZ-minZ << std::endl;
610 // uint16_t *img = new uint16_t[int(maxX+0.5)*int(maxY+0.5)];
611 uint16_t *img = new uint16_t[int(maxX*4.)*int(maxY*4.)];
613 // Set whole image to 0
614 for(int i3=0;i3<int(maxX*4.)*int(maxY*4.);i3++)
615 // for(int i3=0;i3<int(maxX+0.5)*int(maxY+0.5);i3++)
618 for(int i2=0; i2<NP; i2++) {
620 int ordX = int(X[i2]*4.-30);
621 int ordY = int(maxY*4.) - int(Y[i2]*4.)+30;
622 // 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);
624 // img[ /*int(maxX) -*/ int(X[i2]*4.-30) + (int(maxY*4.) - int(Y[i2]*4.)+30) * int(maxX*4.) ] = int(tabVal[i2]*100);
625 img[ /*int(maxX) -*/ ordX + ordY * int(maxX*4.) ] = int(tabVal[i2]*100);
627 // Try to round up, just to see.
628 for(int iii=ordY-3; iii<ordY+4; iii++) for(int jjj=ordX-3; jjj<ordX+4; jjj++)
629 img[ jjj + iii * int(maxX*4.) ] = int(tabVal[i2]*100);
630 std::cout << int(X[i2]*4.) << " " << int(Y[i2]*4.) << " = " << int(tabVal[i2]*100) << std::endl;
633 // GDCM_NAME_SPACE::Debug::DebugOn();
635 std::ostringstream str;
637 GDCM_NAME_SPACE::File *file;
638 file = GDCM_NAME_SPACE::File::New();
640 // Set the image size
642 str << (int)(maxX*4.);
643 file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
645 str << (int)(maxY*4.);
646 file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
648 // Set the pixel type
650 file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
652 str << 16; // may be 12 or 16 if componentSize =16
653 file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
654 file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
656 // Set the pixel representation // 0/1 , 0=unsigned
657 file->InsertEntryString("1",0x0028,0x0103, "US"); // Pixel Representation
658 // Set the samples per pixel // 1:Grey level, 3:RGB
659 file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
662 // Set Rescale Intercept
665 file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
670 file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
673 GDCM_NAME_SPACE::FileHelper *fileH;
674 fileH = GDCM_NAME_SPACE::FileHelper::New(file);
675 // cast is just to avoid warnings (*no* conversion)
676 //fileH->SetImageData((uint8_t *)img,int(maxX*maxY)*sizeof(uint16_t)); // troubles when maxX, mayY are *actually* float!
677 fileH->SetImageData((uint8_t *)img,int(maxX*4.)*int(maxY*4.)*sizeof(uint16_t));
678 fileH->SetWriteModeToRaw();
679 fileH->SetWriteTypeToDcmExplVR();
681 if( !fileH->Write(dcmImageName))
682 std::cout << "Failed for [" << dcmImageName << "]\n"
683 << " File is unwrittable" << std::endl;