[Ensure slice timing always encoded for Siemens EPI](neurolabusc/dcm_…

…qa#4 (comment))

COMPILE.md

@@ -140,4 +140,4 @@ CXXFLAGS = -DmyUseCOut -DWX_PRECOMP ....
  i.) For a full refresh
 rm clipboard
 rm *.o
-make
+make

README.md

@@ -11,6 +11,9 @@ This software is open source. The bulk of the code is covered by the BSD license
 
 ## Versions
 
+23-June-2017
+ - [Ensure slice timing always encoded for Siemens EPI](https://github.com/neurolabusc/dcm_qa/issues/4#issuecomment-310707906)
+
 21-June-2017
  - Read DICOM header in 1Mb segments rather than loading whole file : reduces ram usage and [faster for systems with slow io](https://github.com/rordenlab/dcm2niix/issues/104).
  - Report [TotalReadoutTime](https://github.com/rordenlab/dcm2niix/issues/98).

console/nii_dicom.cpp

@@ -621,6 +621,7 @@ struct TDICOMdata clear_dicom_data() {
         d.angulation[i] = 0.0f;
         d.xyzMM[i] = 1;
     }
+    d.CSA.sliceTiming[0] = -1.0f; //impossible value denotes not known
     d.CSA.numDti = 0;
     for (int i=0; i < 5; i++)
         d.xyzDim[i] = 1;
@@ -1004,19 +1005,16 @@ int readCSAImageHeader(unsigned char *buff, int lLength, struct TCSAdata *CSA, i
             else if (strcmp(tagCSA.name, "BandwidthPerPixelPhaseEncode") == 0)
                 CSA->bandwidthPerPixelPhaseEncode = csaMultiFloat (&buff[lPos], 3,lFloats, &itemsOK);
             else if ((strcmp(tagCSA.name, "MosaicRefAcqTimes") == 0) && (tagCSA.nitems > 3)  ){
-//#ifdef _MSC_VER
 				float * sliceTimes = (float *)malloc(sizeof(float) * (tagCSA.nitems + 1));
-//#else
-//				float sliceTimes[tagCSA.nitems + 1];
-//#endif
                 csaMultiFloat (&buff[lPos], tagCSA.nitems,sliceTimes, &itemsOK);
                 float maxTimeValue, minTimeValue, timeValue1;
-                for (int z = 0; z < kMaxDTI4D; z++)
-        			dti4D->S[z].sliceTiming = -1.0;
-
-                if (itemsOK <= kMaxDTI4D)
+                for (int z = 0; z < kMaxEPI3D; z++)
+        			CSA->sliceTiming[z] = -1.0;
+        		if (itemsOK <= kMaxEPI3D) {
                 	for (int z = 1; z <= itemsOK; z++)
-                		dti4D->S[z-1].sliceTiming = sliceTimes[z];
+                		CSA->sliceTiming[z-1] = sliceTimes[z];
+                } else
+                	printError("Please increase kMaxEPI3D and recompile\n");
                 CSA->multiBandFactor = 1;
                 timeValue1 = sliceTimes[1];
                 int nTimeZero = 0;
@@ -1328,10 +1326,8 @@ struct TDICOMdata  nii_readParRec (char * parname, int isVerbose, struct TDTI4D
     char *p = fgets (buff, LINESZ, fp);
     bool isIntenScaleVaries = false;
     bool isIndexSequential = true;
-    for (int i = 0; i < kMaxDTI4D; i++) {
+    for (int i = 0; i < kMaxDTI4D; i++)
         dti4D->S[i].V[0] = -1.0;
-        dti4D->S[i].sliceTiming = -1.0;
-    }
     //d.dti4D = (TDTI *)malloc(kMaxDTI4D * sizeof(TDTI));
     while (p) {
         if (strlen(buff) < 1)
@@ -2650,7 +2646,6 @@ struct TDICOMdata readDICOMv(char * fname, int isVerbose, int compressFlag, stru
 #define kNest 0xFFFE +(0xE000 << 16 ) //Item follows SQ
 #define  kUnnest 0xFFFE +(0xE00D << 16 ) //ItemDelimitationItem [length defined] http://www.dabsoft.ch/dicom/5/7.5/
 #define  kUnnest2 0xFFFE +(0xE0DD << 16 )//SequenceDelimitationItem [length undefined]
-    dti4D->S[0].sliceTiming = -1.0;
     int nest = 0;
     double zSpacing = -1.0l; //includes slice thickness plus gap
     int locationsInAcquisitionGE = 0;

console/nii_dicom.h

@@ -38,8 +38,9 @@ extern "C" {
 	#define kCCsuf " CompilerNA" //unknown compiler!
 #endif
 
- #define kDCMvers "v1.0.20170621" kDCMsuf kCCsuf
+ #define kDCMvers "v1.0.20170623" kDCMsuf kCCsuf
 
+static const int kMaxEPI3D = 1024; //maximum number of EPI images in Siemens Mosaic
 static const int kMaxDTI4D = 4096; //maximum number of DTI directions for 4D (Philips) images, also maximum number of 3D slices for Philips 3D and 4D images
 #define kDICOMStr 64
 #define kMANUFACTURER_UNKNOWN  0
@@ -58,7 +59,6 @@ static const int kCompressC3 = 2; //obsolete JPEG lossless
 static const int kCompress50 = 3; //obsolete JPEG lossy
     struct TDTI {
         float V[4];
-        float sliceTiming;
         int sliceNumberMrPhilips;
     };
     struct TDTI4D {
@@ -93,7 +93,7 @@ static const int kCompress50 = 3; //obsolete JPEG lossy
 #endif
     struct TCSAdata {
     	bool isPhaseMap;
-        float dtiV[4], sliceNormV[4], bandwidthPerPixelPhaseEncode, sliceMeasurementDuration;
+        float sliceTiming[kMaxEPI3D], dtiV[4], sliceNormV[4], bandwidthPerPixelPhaseEncode, sliceMeasurementDuration;
         int numDti, SeriesHeader_offset, SeriesHeader_length, multiBandFactor, sliceOrder, slice_start, slice_end, mosaicSlices,protocolSliceNumber1,phaseEncodingDirectionPositive;
     };
     struct TDICOMdata {

console/nii_dicom_batch.cpp

@@ -664,17 +664,13 @@ void nii_SaveBIDS(char pathoutname[], struct TDICOMdata d, struct TDCMopts opts,
     		if (effectiveEchoSpacing > 0.0)
     			fprintf(fp, "\t\"TrueEchoSpacing\": %g,\n", effectiveEchoSpacing * d.accelFactPE);
 	}
-	bool first = 1;
-	if (dti4D->S[0].sliceTiming >= 0.0) {
+	if (d.CSA.sliceTiming[0] >= 0.0) {
    		fprintf(fp, "\t\"SliceTiming\": [\n");
-   		for (int i = 0; i < kMaxDTI4D; i++) {
-			if (dti4D->S[i].sliceTiming >= 0.0){
-			  if (!first)
-				  fprintf(fp, ",\n");
-				else
-				  first = 0;
-				fprintf(fp, "\t\t%g", dti4D->S[i].sliceTiming / 1000.0 );
-			}
+   		for (int i = 0; i < kMaxEPI3D; i++) {
+   			if (d.CSA.sliceTiming[i] < 0.0) break;
+			if (i != 0)
+				fprintf(fp, ",\n");
+			fprintf(fp, "\t\t%g", d.CSA.sliceTiming[i] / 1000.0 );
 		}
 		fprintf(fp, "\t],\n");
 	}

console/nii_foreign.cpp

@@ -400,7 +400,6 @@ int  convert_foreign (const char *fn, struct TDCMopts opts){
 	printMessage("Saving ECAT as '%s'\n", niiFilename);
 	if (ret != EXIT_SUCCESS) return ret;
 	struct TDTI4D dti4D;
-	dti4D.S[0].sliceTiming = -1.0;
 	nii_SaveBIDS(niiFilename, dcm, opts, &dti4D, &hdr, fn);
 	ret = nii_saveNII(niiFilename, hdr, img, opts);
 	free(img);