Skip to content

Commit

Permalink
Experimental BEP001 support (bids-standard/bids-specification#681)
Browse files Browse the repository at this point in the history
  • Loading branch information
@neurolabusc
neurolabusc committed Jun 16, 2021
1 parent 4fc19d0 commit 9af237f
Show file tree
Hide file tree
Showing 3 changed files with 47 additions and 24 deletions.
14 changes: 13 additions & 1 deletion console/nii_dicom.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -853,6 +853,7 @@ struct TDICOMdata clear_dicom_data() {
d.durationLabelPulseGE = -1;
d.aslFlagsGE = 0;
d.partialFourierDirection = kPARTIAL_FOURIER_DIRECTION_UNKNOWN;
d.mtState = -1;
d.interp3D = -1;
for (int i = 0; i < kMaxOverlay; i++)
d.overlayStart[i] = 0;
Expand Down Expand Up @@ -4257,6 +4258,8 @@ struct TDICOMdata readDICOMx(char * fname, struct TDCMprefs* prefs, struct TDTI4
#define kPatientOrient 0x0018+(0x5100<< 16 ) //0018,5100. patient orientation - 'HFS'
#define kInversionRecovery 0x0018+uint32_t(0x9009 << 16 ) //'CS' 'YES'/'NO'
#define kEchoPlanarPulseSequence 0x0018+uint32_t(0x9018 << 16 ) //'CS' 'YES'/'NO'
#define kMagnetizationTransferAttribute 0x0018+uint32_t(0x9020 << 16 ) //'CS' 'ON_RESONANCE','OFF_RESONANCE','NONE'

#define kRectilinearPhaseEncodeReordering 0x0018+uint32_t(0x9034 << 16) //'CS' 'REVERSE_LINEAR'/'LINEAR'
#define kPartialFourierDirection 0x0018+uint32_t(0x9036 << 16) //'CS'
#define kParallelReductionFactorInPlane 0x0018+uint32_t(0x9069<< 16 ) //FD
Expand Down Expand Up @@ -5379,7 +5382,16 @@ uint32_t kSequenceDelimitationItemTag = 0xFFFE +(0xE0DD << 16 );
if (toupper(buffer[lPos]) == 'Y')
d.isEPI = true;
break;
case kRectilinearPhaseEncodeReordering : { //'CS' [REVERSE_LINEAR],[LINEAR],[CENTRIC],[REVERSE_CENTRIC]
case kMagnetizationTransferAttribute : //'CS' 'ON_RESONANCE','OFF_RESONANCE','NONE'
if (lLength < 2) break; //https://github.com/bids-standard/bids-specification/pull/681
if ((toupper(buffer[lPos]) == 'O') && (toupper(buffer[lPos]) == 'F')) // OFF_RESONANCE
d.mtState = 1; //TRUE
if ((toupper(buffer[lPos]) == 'O') && (toupper(buffer[lPos]) == 'N')) // ON_RESONANCE and NONE
d.mtState = 0; //FALSE
if ((toupper(buffer[lPos]) == 'N') && (toupper(buffer[lPos]) == 'O')) // ON_RESONANCE and NONE
d.mtState = 0; //FALSE
break;
case kRectilinearPhaseEncodeReordering : { //'CS' [REVERSE_LINEAR],[LINEAR],[CENTRIC],[REVERSE_CENTRIC]
if (d.manufacturer != kMANUFACTURER_GE) break; //only found in GE software beginning with RX27
if (lLength < 2) break;
if (toupper(buffer[lPos]) == 'L')
Expand Down
4 changes: 2 additions & 2 deletions console/nii_dicom.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ extern "C" {
#define kCPUsuf " " //unknown CPU
#endif

#define kDCMdate "v1.0.20210531"
#define kDCMdate "v1.0.20210616"
#define kDCMvers kDCMdate " " kJP2suf kLSsuf kCCsuf kCPUsuf

static const int kMaxEPI3D = 1024; //maximum number of EPI images in Siemens Mosaic
Expand Down Expand Up @@ -194,7 +194,7 @@ static const uint8_t MAX_NUMBER_OF_DIMENSIONS = 8;
int xyzDim[5];
uint32_t coilCrc, seriesUidCrc, instanceUidCrc;
int overlayStart[kMaxOverlay];
int partialFourierDirection, interp3D, aslFlagsGE, durationLabelPulseGE, epiVersionGE, internalepiVersionGE, maxEchoNumGE, rawDataRunNumber, numberOfImagesInGridUIH, numberOfDiffusionDirectionGE, phaseEncodingGE, protocolBlockStartGE, protocolBlockLengthGE, modality, dwellTime, effectiveEchoSpacingGE, phaseEncodingLines, phaseEncodingSteps, echoTrainLength, echoNum, sliceOrient, manufacturer, converted2NII, acquNum, imageNum, imageStart, imageBytes, bitsStored, bitsAllocated, samplesPerPixel,locationsInAcquisition, locationsInAcquisitionConflict, compressionScheme;
int mtState, partialFourierDirection, interp3D, aslFlagsGE, durationLabelPulseGE, epiVersionGE, internalepiVersionGE, maxEchoNumGE, rawDataRunNumber, numberOfImagesInGridUIH, numberOfDiffusionDirectionGE, phaseEncodingGE, protocolBlockStartGE, protocolBlockLengthGE, modality, dwellTime, effectiveEchoSpacingGE, phaseEncodingLines, phaseEncodingSteps, echoTrainLength, echoNum, sliceOrient, manufacturer, converted2NII, acquNum, imageNum, imageStart, imageBytes, bitsStored, bitsAllocated, samplesPerPixel,locationsInAcquisition, locationsInAcquisitionConflict, compressionScheme;
float groupDelay, decayFactor, percentSampling,waterFatShift, numberOfAverages, imagingFrequency, patientWeight, zSpacing, zThick, pixelBandwidth, SAR, phaseFieldofView, accelFactPE, accelFactOOP, flipAngle, fieldStrength, TE, TI, TR, intenScale, intenIntercept, intenScalePhilips, gantryTilt, lastScanLoc, angulation[4];
float orient[7], patientPosition[4], patientPositionLast[4], xyzMM[4], stackOffcentre[4];
float rtia_timerGE, radionuclidePositronFraction, radionuclideTotalDose, radionuclideHalfLife, doseCalibrationFactor; //PET ISOTOPE MODULE ATTRIBUTES (C.8-57)
Expand Down
53 changes: 32 additions & 21 deletions console/nii_dicom_batch.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -570,7 +570,7 @@ int phoenixOffsetCSASeriesHeader(unsigned char *buff, int lLength) {
typedef struct {
float TE0, TE1, delayTimeInTR, phaseOversampling, phaseResolution, txRefAmp;
int phaseEncodingLines, existUcImageNumb, ucMode, baseResolution, interp, partialFourier,echoSpacing,
difBipolar, parallelReductionFactorInPlane, refLinesPE;
difBipolar, parallelReductionFactorInPlane, refLinesPE, combineMode, patMode;
float alFree[kMaxWipFree] ;
float adFree[kMaxWipFree];
float alTI[kMaxWipFree];
Expand All @@ -596,6 +596,8 @@ void siemensCsaAscii(const char * filename, TCsaAscii* csaAscii, int csaOffset,
csaAscii->difBipolar = 0; //0=not assigned,1=bipolar,2=monopolar
csaAscii->parallelReductionFactorInPlane = 0;
csaAscii->refLinesPE = 0;
csaAscii->combineMode = 0;
csaAscii->patMode = 0;
for (int i = 0; i < 8; i++)
shimSetting[i] = 0.0;
strcpy(coilID, "");
Expand Down Expand Up @@ -646,22 +648,6 @@ void siemensCsaAscii(const char * filename, TCsaAscii* csaAscii, int csaOffset,
csaAscii->phaseEncodingLines = readKey(keyStrLns, keyPos, csaLengthTrim);
char keyStrUcImg[] = "sSliceArray.ucImageNumb";
csaAscii->existUcImageNumb = readKey(keyStrUcImg, keyPos, csaLengthTrim);
/*
//TODO!
int combineMode = -1;
char keyStrCombineMode[] = "ucCoilCombineMode";
combineMode = readKeyN1(keyStrCombineMode, keyPos, csaLengthTrim);
//BIDS CoilCombinationMethod <- Siemens 'Coil Combine Mode' CSA ucCoilCombineMode 1 = Sum of Squares, 2 = Adaptive Combine,
printf("CoilCombineMode %d\n", combineMode);
//BIDS MatrixCoilMode <- Siemens 'PAT mode' CSA ucPATMode ?? 1?= "SENSE" 2= "GRAPPA",
//sPat.ucPATMode = 2
int patMode = -1;
char keyStrPATMode[] = "sPat.ucPATMode"; //n.b. field set even if PAT not enabled, e.g. will list SENSE for a R-factor of 1
patMode = readKeyN1(keyStrPATMode, keyPos, csaLengthTrim);
printf("PATMODE %d\n", patMode);
*/
char keyStrUcMode[] = "sSliceArray.ucMode";
csaAscii->ucMode = readKeyN1(keyStrUcMode, keyPos, csaLengthTrim);
char keyStrBase[] = "sKSpace.lBaseResolution";
Expand All @@ -685,6 +671,13 @@ void siemensCsaAscii(const char * filename, TCsaAscii* csaAscii, int csaOffset,
csaAscii->parallelReductionFactorInPlane = readKey(keyStrAF, keyPos, csaLengthTrim);
char keyStrRef[] = "sPat.lRefLinesPE";
csaAscii->refLinesPE = readKey(keyStrRef, keyPos, csaLengthTrim);
char keyStrCombineMode[] = "ucCoilCombineMode";
csaAscii->combineMode = readKeyN1(keyStrCombineMode, keyPos, csaLengthTrim);
//BIDS CoilCombinationMethod <- Siemens 'Coil Combine Mode' CSA ucCoilCombineMode 1 = Sum of Squares, 2 = Adaptive Combine,
//printf("CoilCombineMode %d\n", csaAscii->combineMode);
char keyStrPATMode[] = "sPat.ucPATMode"; //n.b. field set even if PAT not enabled, e.g. will list SENSE for a R-factor of 1
csaAscii->patMode = readKeyN1(keyStrPATMode, keyPos, csaLengthTrim);
//printf("PATMODE %d\n", csaAscii->patMode);
//char keyStrETD[] = "sFastImaging.lEchoTrainDuration";
//*echoTrainDuration = readKey(keyStrETD, keyPos, csaLengthTrim);
//char keyStrEF[] = "sFastImaging.lEPIFactor";
Expand Down Expand Up @@ -984,6 +977,17 @@ void json_Float(FILE *fp, const char *sLabel, float sVal) {
fprintf(fp, sLabel, sVal );
} //json_Float

void json_Bool(FILE *fp, const char *sLabel, int sVal) {
// json_Str(fp, "\t\"MTState\"", d.mtState);
//n.b. in JSON, true and false are lower case, whereas in Python they are capitalized
// only print 0 / 1 for false true, ignore negative values
if (sVal == 0) fprintf(fp, sLabel, "false");
if (sVal == 1) fprintf(fp, sLabel, "true");

//if (sVal = 0) fprintf(" : false,\n" );
//if (sVal = 1) fprintf(" : true,\n" );
} //json_Bool

void rescueProtocolName(struct TDICOMdata *d, const char * filename) {
//tools like gdcmanon strip protocol name (0018,1030) but for Siemens we can recover it from CSASeriesHeaderInfo (0029,1020)
if ((d->manufacturer != kMANUFACTURER_SIEMENS) || (d->CSA.SeriesHeader_offset < 1) || (d->CSA.SeriesHeader_length < 1)) return;
Expand Down Expand Up @@ -1263,9 +1267,7 @@ tse3d: T2*/
fprintf(fp, "\t\t%g", dti4D->frameDuration[i] / 1000.0 ); // from 0018,1242 ms -> sec
}
fprintf(fp, "\t],\n");
}


}
//CT parameters
if ((d.TE > 0.0) && (d.isXRay)) fprintf(fp, "\t\"XRayExposure\": %g,\n", d.TE );
//MRI parameters
Expand All @@ -1282,6 +1284,7 @@ tse3d: T2*/
json_Float(fp, "\t\"RepetitionTime\": %g,\n", d.TR / 1000.0 );
json_Float(fp, "\t\"RepetitionTimeExcitation\": %g,\n", dti4D->repetitionTimeExcitation);
json_Float(fp, "\t\"RepetitionTimeInversion\": %g,\n", dti4D->repetitionTimeInversion);
json_Bool(fp, "\t\"MTState\": %s,\n", d.mtState);
json_Float(fp, "\t\"InversionTime\": %g,\n", d.TI / 1000.0 );
json_Float(fp, "\t\"FlipAngle\": %g,\n", d.flipAngle );
bool interp = false; //2D interpolation
Expand Down Expand Up @@ -1502,8 +1505,17 @@ tse3d: T2*/
json_Str(fp, "\t\"ProtocolName\": \"%s\",\n", protocolName);
if (csaAscii.refLinesPE > 0)
fprintf(fp, "\t\"RefLinesPE\": %d,\n", csaAscii.refLinesPE);
//https://github.com/bids-standard/bids-specification/pull/681#issuecomment-861767213
if (csaAscii.combineMode == 1)
fprintf(fp, "\t\"CoilCombinationMethod\": \"Sum of Squares\",\n" );
if (csaAscii.combineMode == 2)
fprintf(fp, "\t\"CoilCombinationMethod\": \"Adaptive Combine\",\n" );
json_Str(fp, "\t\"ConsistencyInfo\": \"%s\",\n", consistencyInfo);
if (csaAscii.parallelReductionFactorInPlane > 0) {//AccelFactorPE -> phase encoding
if (csaAscii.patMode == 1)
fprintf(fp, "\t\"MatrixCoilMode\": \"SENSE\",\n" );
if (csaAscii.patMode == 2)
fprintf(fp, "\t\"MatrixCoilMode\": \"GRAPPA\",\n" );
if (d.accelFactPE < 1.0) { //value not found in DICOM header, but WAS found in CSA ascii
d.accelFactPE = csaAscii.parallelReductionFactorInPlane; //value found in ASCII but not in DICOM (0051,1011)
//fprintf(fp, "\t\"ParallelReductionFactorInPlane\": %g,\n", d.accelFactPE);
Expand All @@ -1529,7 +1541,6 @@ tse3d: T2*/
}
#endif
//GE ASL specific tags

if (d.aslFlagsGE & kASL_FLAG_GE_CONTINUOUS)
fprintf(fp, "\t\"ASLContrastTechnique\": \"CONTINUOUS\",\n" );
if (d.aslFlagsGE & kASL_FLAG_GE_PSEUDOCONTINUOUS)
Expand Down

0 comments on commit 9af237f

Please sign in to comment.