PET image index and variable timing #184

rordenlab · Apr 26, 2018 · b699b2b · b699b2b
1 parent c3188e3
commit b699b2b
Show file tree

Hide file tree

Showing 2 changed files with 38 additions and 0 deletions.
diff --git a/console/nii_dicom.cpp b/console/nii_dicom.cpp
@@ -3403,6 +3403,7 @@ const uint32_t kEffectiveTE  = 0x0018+ (0x9082 << 16);
 //If ImageType is REPROJECTION we slice direction is reversed - need example to test
 // #define  kSeriesType  0x0054+(0x1000 << 16 )
 #define  kDoseCalibrationFactor  0x0054+(0x1322<< 16 )
+#define  kPETImageIndex  0x0054+(0x1330<< 16 )
 #define  kIconImageSequence 0x0088+(0x0200 << 16 )
 #define  kDiffusionBFactor  0x2001+(0x1003 << 16 )// FL
 #define  kSliceNumberMrPhilips 0x2001+(0x100A << 16 ) //IS Slice_Number_MR
@@ -3434,6 +3435,7 @@ uint32_t kSequenceDelimitationItemTag = 0xFFFE +(0xE0DD << 16 );
 double TE = 0.0; //most recent echo time recorded
 	bool is2005140FSQ = false;
     int locationsInAcquisitionGE = 0;
+    int PETImageIndex = 0;
     int inStackPositionNumber = 0;
     int maxInStackPositionNumber = 0;
     int temporalPositionIdentifier = 0;
@@ -4087,6 +4089,9 @@ double TE = 0.0; //most recent echo time recorded
             case kDoseCalibrationFactor :
                 d.doseCalibrationFactor = dcmStrFloat(lLength, &buffer[lPos]);
                 break;
+            case kPETImageIndex :
+            	PETImageIndex = dcmInt(lLength,&buffer[lPos],d.isLittleEndian);
+            	break;
             case kBitsAllocated :
                 d.bitsAllocated = dcmInt(lLength,&buffer[lPos],d.isLittleEndian);
                 break;
@@ -4622,6 +4627,10 @@ double TE = 0.0; //most recent echo time recorded
     	printWarning("DICOM violation (contact vendor): compressed image without image fragments, assuming image offset defined by 0x7FE0,x0010: %s\n", fname);
     	d.imageStart = encapsulatedDataImageStart;
     }
+    if ((d.modality == kMODALITY_PT) && (PETImageIndex > 0)) {
+    	d.imageNum = PETImageIndex; //https://github.com/rordenlab/dcm2niix/issues/184
+    	//printWarning("PET scan using 0054,1330 for image number %d\n", PETImageIndex);
+    }
     //Recent Philips images include DateTime (0008,002A) but not separate date and time (0008,0022 and 0008,0032)
     #define kYYYYMMDDlen 8 //how many characters to encode year,month,day in "YYYYDDMM" format
     if ((strlen(acquisitionDateTimeTxt) > (kYYYYMMDDlen+5)) && (!isFloatDiff(d.acquisitionTime, 0.0f)) && (!isFloatDiff(d.acquisitionDate, 0.0f)) ) {

diff --git a/console/nii_dicom_batch.cpp b/console/nii_dicom_batch.cpp
@@ -2466,6 +2466,10 @@ int nii_saveCrop(char * niiFilename, struct nifti_1_header hdr, unsigned char* i
     return returnCode;
 }// nii_saveCrop()
 
+float acquisitionTimeDifference(struct TDICOMdata * d, struct TDICOMdata * d1) {
+	if (d->acquisitionDate != d1->acquisitionDate) return -1; //to do: scans running across midnight
+	return (d1->acquisitionTime - d->acquisitionTime);
+}
 void checkDateTimeOrder(struct TDICOMdata * d, struct TDICOMdata * d1) {
 	if (d->acquisitionDate < d1->acquisitionDate) return; //d1 occurred on later date
 	if (d->acquisitionTime <= d1->acquisitionTime) return; //d1 occurred on later (or same) time
@@ -2579,6 +2583,31 @@ int saveDcm2NiiCore(int nConvert, struct TDCMsort dcmSort[],struct TDICOMdata dc
                     printMessage("Slice positions repeated, but number of slices (%d) not divisible by number of repeats (%d): missing images?\n", nConvert, nAcq);
                 }
             }
+            //next: detect variable inter-volume time https://github.com/rordenlab/dcm2niix/issues/184
+    		bool trVaries = false;
+    		float tr = -1;
+    		int prevVolIndx = indx0;
+    		for (int i = 0; i < nConvert; i++)
+                    if (isSamePosition(dcmList[indx0],dcmList[dcmSort[i].indx])) {
+                    	float trDiff = acquisitionTimeDifference(&dcmList[prevVolIndx], &dcmList[dcmSort[i].indx]);
+                    	prevVolIndx = dcmSort[i].indx;
+                    	if (trDiff <= 0) continue;
+                    	if (tr < 0) tr = trDiff;
+                    	if (!isSameFloatGE(tr,trDiff))
+                    		trVaries = true;
+            		}
+            if (trVaries) {
+            	saveAs3D = true;
+            	printWarning("Creating independent volumes as time between volumes varies\n");
+				printMessage(" OnsetTime = [");
+				for (int i = 0; i < nConvert; i++)
+						if (isSamePosition(dcmList[indx0],dcmList[dcmSort[i].indx])) {
+							float trDiff = acquisitionTimeDifference(&dcmList[indx0], &dcmList[dcmSort[i].indx]);
+							printMessage(" %g", trDiff);
+						}
+				printMessage(" ]\n");
+            }
+            //next: detect variable inter-slice distance
             float dx = intersliceDistance(dcmList[dcmSort[0].indx],dcmList[dcmSort[1].indx]);
             bool dxVaries = false;
             for (int i = 1; i < nConvert; i++)