Merge remote branch 'origin/master'

ANTsX · Apr 1, 2017 · 646459b · 646459b
2 parents 1df1128 + a26cf76
commit 646459b
Show file tree

Hide file tree

Showing 5 changed files with 154 additions and 99 deletions.
diff --git a/Examples/antsApplyTransforms.cxx b/Examples/antsApplyTransforms.cxx
@@ -109,6 +109,31 @@ unsigned int numTensorElements<0>()
   return 0;
 }
 
+template <unsigned int NDim>
+std::vector<unsigned int> tensorDiagonalArrayIndices()
+{
+  std::vector<unsigned int> diagElements;
+  for( unsigned int d = 0; d < NDim; d++ )
+    { // I think this is correct for upper-triangular ordering but only tested on 3D tensors
+    diagElements.push_back(d * NDim - d * (d - 1) / 2);
+    }
+  return diagElements;
+}
+
+template <unsigned int NDim>
+bool isDiagonalElement(std::vector<unsigned int> diagElements, unsigned int ind)
+{
+  for( unsigned int i = 0; i < NDim; i++ )
+    {
+    if ( diagElements[i] == ind )
+      {
+      return true;
+      }
+    }
+
+  return false;
+}
+
 template <class T, unsigned int Dimension>
 int antsApplyTransforms( itk::ants::CommandLineParser::Pointer & parser, unsigned int inputImageType = 0 )
 {
@@ -133,6 +158,8 @@ int antsApplyTransforms( itk::ants::CommandLineParser::Pointer & parser, unsigne
 
   const unsigned int NumberOfTensorElements = numTensorElements<Dimension>();
 
+  std::vector<unsigned int> tensorDiagIndices = tensorDiagonalArrayIndices<Dimension>();
+
   typename TimeSeriesImageType::Pointer timeSeriesImage = ITK_NULLPTR;
   typename TensorImageType::Pointer tensorImage = ITK_NULLPTR;
   typename DisplacementFieldType::Pointer vectorImage = ITK_NULLPTR;
@@ -379,7 +406,10 @@ int antsApplyTransforms( itk::ants::CommandLineParser::Pointer & parser, unsigne
     }
   if( verbose )
     {
-    std::cout << "Default pixel value: " << defaultValue << std::endl;
+    if( inputImageType == 2 )
+      std::cout << "Default pixel mean diffusivity: " << defaultValue << std::endl;
+    else
+      std::cout << "Default pixel value: " << defaultValue << std::endl;
     }
   for( unsigned int n = 0; n < inputImages.size(); n++ )
     {
@@ -388,7 +418,24 @@ int antsApplyTransforms( itk::ants::CommandLineParser::Pointer & parser, unsigne
     resampleFilter->SetInput( inputImages[n] );
     resampleFilter->SetOutputParametersFromImage( referenceImage );
     resampleFilter->SetTransform( compositeTransform );
-    resampleFilter->SetDefaultPixelValue( defaultValue );
+    if ( inputImageType == 2 )
+      {
+      // Set background pixel values in tensor images to produce an isotropic tensor
+      if ( defaultValue > 0 && isDiagonalElement<Dimension>(tensorDiagIndices, n) )
+        {
+        // defaultValue == MD of isotropic tensor. Resampling is done in log space
+        resampleFilter->SetDefaultPixelValue( log( defaultValue ) );
+        }
+             else
+        {
+        resampleFilter->SetDefaultPixelValue( 0 );
+        }
+      }
+    else
+      {
+      // for non-tensor images, set the same background value for each component
+      resampleFilter->SetDefaultPixelValue( defaultValue );
+      }
 
     interpolator->SetInputImage( inputImages[n] );
     resampleFilter->SetInterpolator( interpolator );
@@ -790,7 +837,8 @@ static void antsApplyTransformsInitializeCommandLineOptions( itk::ants::CommandL
   std::string description =
     std::string( "Default voxel value to be used with input images only. " )
     + std::string( "Specifies the voxel value when the input point maps outside " )
-    + std::string( "the output domain" );
+    + std::string( "the output domain. With tensor input images, specifies the " )
+    + std::string( "default voxel eigenvalues. " );
 
   OptionType::Pointer option = OptionType::New();
   option->SetLongName( "default-value" );

diff --git a/Scripts/antsASLProcessing.R b/Scripts/antsASLProcessing.R
@@ -21,15 +21,15 @@ optlist <- list(
               help='Full width half max for smoothing'),
   make_option(c('-m', '--method'), default='regression',
               help=paste(' method for perfusion calculation. \n\t\tOne of:',
-                '"SimpleSubtraction", "SurroundSubtraction", "SincSubtraction",',
+                '"regression", "subtraction", "bayesian",',
                 '"RobustRegression", "BayesianRegression", "LocalBayesianRegression."')),
   make_option(c('-d', '--denoising'), default='CompCorMotion',
                 help=paste('denoising method.',
                 'Options are: "CompCor", "Motion", "Detrending",',
                 '\n\t\t"Cross-Validation", "OutlierRejection".',
                 'Multiple options can be specified',
                 '(e.g., "CompCorMotion" is legal).  Default is %default.')),
-  make_option(c('-g', '--debug'), default=FALSE, action='store_true',
+  make_option(c('-g', '--debug'), default=0,
                help=paste('Save debugging information, including motion',
                   'correction and nuisance variables')),
   make_option(c('-b', '--bloodT1'), default=0.67,
@@ -40,13 +40,13 @@ optlist <- list(
               help=' number of bootstrap samples (defaults to %default)'),
   make_option(c('-e', '--bootstrapPercent'), default=0.70,
               help='percent to sample per bootstrap run (defaults to %default)'),
-  make_option(c('-k', '--keepTmp'), default=F, action='store_true',
+  make_option(c('-k', '--keepTmp'), default=0,
               help=paste('keep tmp files, including warps',
                          '(defaults to %default--takes lots of space to save)')),
-  make_option(c('-f', '--bootstrapReplace'), default=F, action='store_true',
+  make_option(c('-f', '--bootstrapReplace'), default=0,
               help=paste('bootstrap with replacement?  takes arguments',
-                         '"false" or "true"; defaults to false.')),
-  make_option(c('-v', '--verbose'), default=F, action='store_true',
+                         '0 or 1; defaults to 0.')),
+  make_option(c('-v', '--verbose'), default=0,
               help='verbose output.'))
 
 usage <- OptionParser(option_list=optlist, usage='Usage: %prog <s> [otlcxmdgbrnekfv]')
@@ -111,7 +111,7 @@ avg <- n3BiasFieldCorrection(avg, 2)
 mask <- getMask(avg, mean(avg), Inf, 2)
 avg[mask==0] <- 0
 
-moco <- antsMotionCalculation(pcasl, fixed=avg, mask=mask,  moreaccurate=4)
+moco <- antsrMotionCalculation(pcasl, fixed=avg, mask=mask)
 tag.first <- config$tagFirst
 ts <- timeseries2matrix(moco$moco_img, moco$moco_mask)
 if (!tag.first) {
@@ -238,7 +238,7 @@ if (nchar(opt$antsCorticalThicknessPrefix) > 0){
                   'does not exist.'))
   })
   reg.t12asl <- antsRegistration(fixed=avg, moving=braint1,
-    typeofTransform="SyNBold", outprefix=as.character(opt$outputpre))
+    typeofTransform="SyNBold" )
   seg.asl <- antsApplyTransforms(avg, seg, reg.t12asl$fwdtransforms, "MultiLabel")
   antsImageWrite(seg.asl, paste(opt$outputpre,
                       "SegmentationWarpedToASL.nii.gz", sep=''))

diff --git a/Scripts/antsLongitudinalCorticalThickness.sh b/Scripts/antsLongitudinalCorticalThickness.sh
@@ -142,6 +142,7 @@ Optional arguments:
      -r:  rigid alignment to SST                This option dictates if the individual subjects are registered to the single
                                                 subject template before running through antsCorticalThickness.  This potentially
                                                 reduces bias caused by subject orientation and voxel spacing (default = 0).
+     -y:  keep temporary files                  Keep brain extraction/segmentation warps, etc (default = 0).
      -z:  Test / debug mode                     If > 0, runs a faster version of the script. Only for testing. Implies -u 0
                                                 in the antsCorticalThickness.sh script (i.e., no random seeding).
                                                 Requires single thread computation for complete reproducibility.
@@ -256,12 +257,13 @@ MALF_LABEL_STRINGS_FOR_PRIORS=()
 ################################################################################
 
 USE_FLOAT_PRECISION=0
+KEEP_TMP_IMAGES=0
 
 if [[ $# -lt 3 ]] ; then
   Usage >&2
   exit 1
 else
-  while getopts "a:b:c:d:e:f:g:h:j:k:l:m:n:o:p:q:r:s:t:u:v:x:w:z:" OPT
+  while getopts "a:b:c:d:e:f:g:h:j:k:l:m:n:o:p:q:r:s:t:u:v:x:w:y:z:" OPT
     do
       case $OPT in
           a)
@@ -341,6 +343,9 @@ else
        ;;
           w) #atropos prior weight for each individual time point
        ATROPOS_SEGMENTATION_PRIOR_WEIGHT_TIMEPOINT=$OPTARG
+       ;;
+          y) # for ants cortical thickness
+       KEEP_TMP_IMAGES=$OPTARG
        ;;
           z) #debug mode
        DEBUG_MODE=$OPTARG
@@ -503,20 +508,42 @@ if [[ ${#ANATOMICAL_IMAGES[@]} -eq ${NUMBER_OF_MODALITIES} ]];
 
     # Won't be quick unless -q 3 was specified
     # But if you are running a longitudinal script without longitudinal data, that may not be the only problem
+    if [[ $DO_REGISTRATION_TO_TEMPLATE -eq 1 ]];
+        then
+    logCmd ${ANTSPATH}/antsCorticalThickness.sh \
+      -d ${DIMENSION}
+      -x ${ATROPOS_SEGMENTATION_INTERNAL_ITERATIONS} \
+      -t ${REGISTRATION_TEMPLATE} \
+      -q ${RUN_FAST_ANTSCT_TO_GROUP_TEMPLATE} \
+      ${SUBJECT_ANATOMICAL_IMAGES} \
+      -e ${BRAIN_TEMPLATE} \
+      -f ${EXTRACTION_REGISTRATION_MASK} \
+      -m ${EXTRACTION_PRIOR} \
+      -k ${KEEP_TMP_IMAGES} \
+      -g ${DENOISE} \
+      -w ${ATROPOS_SEGMENTATION_PRIOR_WEIGHT_TIMEPOINT} \
+      -z ${DEBUG_MODE} \
+      -p ${SEGMENTATION_PRIOR} \
+      -o ${OUTPUT_PREFIX}
+    fi
 
+    if [[ $DO_REGISTRATION_TO_TEMPLATE -eq 0 ]];
+        then
     logCmd ${ANTSPATH}/antsCorticalThickness.sh \
-      -d ${DIMENSION} -x ${ATROPOS_SEGMENTATION_INTERNAL_ITERATIONS} \
+      -d ${DIMENSION}
+      -x ${ATROPOS_SEGMENTATION_INTERNAL_ITERATIONS} \
       -q ${RUN_FAST_ANTSCT_TO_GROUP_TEMPLATE} \
       ${SUBJECT_ANATOMICAL_IMAGES} \
       -e ${BRAIN_TEMPLATE} \
       -f ${EXTRACTION_REGISTRATION_MASK} \
       -m ${EXTRACTION_PRIOR} \
-      -k 0 \
+      -k ${KEEP_TMP_IMAGES} \
       -g ${DENOISE} \
       -w ${ATROPOS_SEGMENTATION_PRIOR_WEIGHT_TIMEPOINT} \
       -z ${DEBUG_MODE} \
       -p ${SEGMENTATION_PRIOR} \
       -o ${OUTPUT_PREFIX}
+    fi
 
     exit 0
   fi
@@ -740,7 +767,7 @@ if [[ ! -f ${SINGLE_SUBJECT_TEMPLATE_CORTICAL_THICKNESS} ]];
           -f ${EXTRACTION_REGISTRATION_MASK} \
           -g ${DENOISE} \
           -m ${EXTRACTION_PRIOR} \
-          -k 0 \
+          -k ${KEEP_TMP_IMAGES} \
           -z ${DEBUG_MODE} \
           -p ${SEGMENTATION_PRIOR} \
           -w ${ATROPOS_SEGMENTATION_PRIOR_WEIGHT_SST} \
@@ -755,7 +782,7 @@ if [[ ! -f ${SINGLE_SUBJECT_TEMPLATE_CORTICAL_THICKNESS} ]];
           -f ${EXTRACTION_REGISTRATION_MASK} \
           -g ${DENOISE} \
           -m ${EXTRACTION_PRIOR} \
-          -k 0 \
+          -k ${KEEP_TMP_IMAGES} \
           -z ${DEBUG_MODE} \
           -p ${SEGMENTATION_PRIOR} \
           -w ${ATROPOS_SEGMENTATION_PRIOR_WEIGHT_SST} \
@@ -973,7 +1000,7 @@ for (( i=0; i < ${#ANATOMICAL_IMAGES[@]}; i+=$NUMBER_OF_MODALITIES ))
       -m ${SINGLE_SUBJECT_TEMPLATE_EXTRACTION_PRIOR} \
       -f ${SINGLE_SUBJECT_TEMPLATE_EXTRACTION_REGISTRATION_MASK} \
       -g ${DENOISE} \
-      -k 0 \
+      -k ${KEEP_TMP_IMAGES} \
       -z ${DEBUG_MODE} \
       -w ${ATROPOS_SEGMENTATION_PRIOR_WEIGHT_TIMEPOINT} \
       -p ${SINGLE_SUBJECT_TEMPLATE_PRIOR} \