Merge pull request #599 from ANTsX/fixRegistrationRandomSeed

Fix registration random seed
ANTsX · Jun 11, 2018 · cf933a0 · cf933a0
2 parents 73cf3b1 + 1c6aaca
commit cf933a0
Show file tree

Hide file tree

Showing 7 changed files with 129 additions and 20 deletions.
diff --git a/Examples/antsAI.cxx b/Examples/antsAI.cxx
@@ -1221,15 +1221,28 @@ int antsAI( itk::ants::CommandLineParser *parser )
     typedef typename itk::Statistics::MersenneTwisterRandomVariateGenerator RandomizerType;
     typename RandomizerType::Pointer randomizer = RandomizerType::New();
 
-    char* antsRandomSeed = getenv( "ANTS_RANDOM_SEED" );
-    if ( antsRandomSeed != NULL )
+    int antsRandomSeed = 1234;
+
+    itk::ants::CommandLineParser::OptionType::Pointer randomSeedOption = parser->GetOption( "random-seed" );
+    if( randomSeedOption && randomSeedOption->GetNumberOfFunctions() )
       {
-      randomizer->SetSeed( atoi( antsRandomSeed ) );
+      antsRandomSeed = parser->Convert<int>( randomSeedOption->GetFunction(0)->GetName() );
       }
     else
       {
-      randomizer->SetSeed( 1234 );
+      char* envSeed = getenv( "ANTS_RANDOM_SEED" );
+
+      if ( envSeed != NULL )
+	{
+	antsRandomSeed = atoi( envSeed );
+	}
+      }
+
+    if ( antsRandomSeed != 0 ) 
+      {
+      randomizer->SetSeed( antsRandomSeed );
       }
+
     unsigned long index = 0;
 
     switch( samplingStrategy )
@@ -1430,7 +1443,7 @@ int antsAI( itk::ants::CommandLineParser *parser )
                   {
                   affineSearchTransform->Scale( bestScale );
                   parametersList.push_back( affineSearchTransform->GetParameters() );
-                }
+		  }
                 else if( strcmp( transform.c_str(), "rigid" ) == 0 )
                   {
                   rigidSearchTransform->SetIdentity();
@@ -1439,19 +1452,19 @@ int antsAI( itk::ants::CommandLineParser *parser )
                   rigidSearchTransform->Translate( searchTranslation, 0 );
                   rigidSearchTransform->SetMatrix( affineSearchTransform->GetMatrix() );
                   parametersList.push_back( rigidSearchTransform->GetParameters() );
-                }
+		  }
                 else if( strcmp( transform.c_str(), "similarity" ) == 0 )
                   {
                   similaritySearchTransform->SetIdentity();
                   similaritySearchTransform->SetCenter( initialTransform->GetCenter() );
                   similaritySearchTransform->SetOffset( initialTransform->GetOffset() );
                   similaritySearchTransform->SetMatrix( affineSearchTransform->GetMatrix() );
                   similaritySearchTransform->SetScale( bestScale );
-
+		  
                   parametersList.push_back( similaritySearchTransform->GetParameters() );
                   }
+		trialCounter++;
                 }
-              trialCounter++;
               }
             }
           }
@@ -1662,6 +1675,19 @@ void InitializeCommandLineOptions( itk::ants::CommandLineParser *parser )
   parser->AddOption( option );
   }
 
+  {
+  std::string description = std::string( "Use a fixed seed for random number generation. " ) 
+    + std::string( "The default fixed seed is overwritten by this value. " )
+    + std::string( "The fixed seed can be any nonzero int value. If the specified seed is zero, " )
+    + std::string( "the system time will be used." );
+
+  OptionType::Pointer option = OptionType::New();
+  option->SetLongName( "random-seed" );
+  option->SetUsageOption( 0, "seedValue" );
+  option->SetDescription( description );
+  parser->AddOption( option );
+  }
+
   {
   std::string description = std::string( "Verbose output." );
 

diff --git a/Examples/antsMotionCorr.cxx b/Examples/antsMotionCorr.cxx
@@ -592,8 +592,25 @@ int ants_motion( itk::ants::CommandLineParser *parser )
       }
     }
 
-  char* antsRandomSeed = getenv( "ANTS_RANDOM_SEED" );
-
+  // Zero seed means use default behavior: registration randomizer seeds from system time
+  // and does not re-seed iterator
+  int antsRandomSeed = 0;
+
+  itk::ants::CommandLineParser::OptionType::Pointer randomSeedOption = parser->GetOption( "random-seed" );
+  if( randomSeedOption && randomSeedOption->GetNumberOfFunctions() )
+    {
+    antsRandomSeed = parser->Convert<int>( randomSeedOption->GetFunction(0)->GetName() );
+    }
+  else
+    {
+    char* envSeed = getenv( "ANTS_RANDOM_SEED" );
+
+    if ( envSeed != NULL )
+      {
+      antsRandomSeed = atoi( envSeed );
+      }
+    }
+
   unsigned int   nparams = 2;
   itk::TimeProbe totalTimer;
   totalTimer.Start();
@@ -1049,9 +1066,9 @@ int ants_motion( itk::ants::CommandLineParser *parser )
       if( std::strcmp( whichTransform.c_str(), "affine" ) == 0 )
         {
         typename AffineRegistrationType::Pointer affineRegistration = AffineRegistrationType::New();
-        if ( antsRandomSeed != NULL )
+        if ( antsRandomSeed != 0 )
           {
-          affineRegistration->MetricSamplingReinitializeSeed( atoi( antsRandomSeed ) );
+          affineRegistration->MetricSamplingReinitializeSeed( antsRandomSeed );
           }
         typename AffineTransformType::Pointer affineTransform = AffineTransformType::New();
         affineTransform->SetIdentity();
@@ -1131,9 +1148,9 @@ int ants_motion( itk::ants::CommandLineParser *parser )
         typedef itk::ImageRegistrationMethodv4<FixedImageType, FixedImageType,
                                                RigidTransformType> RigidRegistrationType;
         typename RigidRegistrationType::Pointer rigidRegistration = RigidRegistrationType::New();
-        if ( antsRandomSeed != NULL )
+        if ( antsRandomSeed != 0 )
           {
-          rigidRegistration->MetricSamplingReinitializeSeed( atoi( antsRandomSeed ) );
+          rigidRegistration->MetricSamplingReinitializeSeed( antsRandomSeed );
           }
         metric->SetFixedImage( preprocessFixedImage );
         metric->SetVirtualDomainFromImage( preprocessFixedImage );
@@ -1800,6 +1817,18 @@ void antsMotionCorrInitializeCommandLineOptions( itk::ants::CommandLineParser *p
   parser->AddOption( option );
   }
 
+  {
+  std::string description = std::string( "Use a fixed seed for random number generation. " ) 
+    + std::string( "By default, the system clock is used to initialize the seeding. " )
+    + std::string( "The fixed seed can be any nonzero int value." );
+  OptionType::Pointer option = OptionType::New();
+  option->SetLongName( "random-seed" );
+  option->SetUsageOption( 0, "seedValue" );
+  option->SetDescription( description );
+  parser->AddOption( option );
+  }
+
+
   {
   std::string description = std::string( "Verbose output." );
 

diff --git a/Examples/antsRegistration.cxx b/Examples/antsRegistration.cxx
@@ -485,6 +485,18 @@ static void antsRegistrationInitializeCommandLineOptions( itk::ants::CommandLine
   parser->AddOption( option );
   }
 
+  {
+  std::string description = std::string( "Use a fixed seed for random number generation. " ) 
+    + std::string( "By default, the system clock is used to initialize the seeding. " )
+    + std::string( "The fixed seed can be any nonzero int value." );
+
+  OptionType::Pointer option = OptionType::New();
+  option->SetLongName( "random-seed" );
+  option->SetUsageOption( 0, "seedValue" );
+  option->SetDescription( description );
+  parser->AddOption( option );
+  }
+
   {
   std::string description = std::string( "Verbose output." );
 

diff --git a/Examples/antsRegistrationTemplateHeader.h b/Examples/antsRegistrationTemplateHeader.h
@@ -58,6 +58,25 @@ DoRegistration(typename ParserType::Pointer & parser)
     regHelper->SetLogStream( cnul );
     }
 
+  OptionType::Pointer fixRandomSeed = parser->GetOption( "random-seed" );
+  if( fixRandomSeed && fixRandomSeed->GetNumberOfFunctions() )
+    {
+    int randomSeed = parser->Convert<int>( fixRandomSeed->GetFunction(0)->GetName() );
+    regHelper->SetRegistrationRandomSeed(randomSeed);
+    }
+  else
+    {
+    char* randomSeedEnv = getenv( "ANTS_RANDOM_SEED" );
+    if ( randomSeedEnv != NULL )
+      {
+      regHelper->SetRegistrationRandomSeed( atoi( randomSeedEnv ) );
+      }
+    else
+      {
+      regHelper->SetRegistrationRandomSeed(0);
+      }
+    }
+
   OptionType::Pointer transformOption = parser->GetOption( "transform" );
   if( !transformOption || transformOption->GetNumberOfFunctions() == 0 )
     {

diff --git a/Examples/itkantsRegistrationHelper.h b/Examples/itkantsRegistrationHelper.h
@@ -699,6 +699,12 @@ class RegistrationHelper : public itk::Object
   itkGetConstMacro( InitializeTransformsPerStage, bool );
   itkBooleanMacro( InitializeTransformsPerStage );
 
+  /**
+   * Set a constant random seed with an int != 0
+   */
+  itkSetMacro( RegistrationRandomSeed, int );
+  itkGetConstMacro( RegistrationRandomSeed, int );
+
   /**
    * turn on winsorize image intensity normalization
    */
@@ -824,11 +830,9 @@ class RegistrationHelper : public itk::Object
     typename RegistrationMethodType::Pointer registrationMethod = RegistrationMethodType::New();
     typedef typename RegistrationMethodType::OutputTransformType  RegistrationMethodTransformType;
 
-    char* antsRandomSeed = getenv( "ANTS_RANDOM_SEED" );
-    if ( antsRandomSeed != NULL )
+    if ( this->m_RegistrationRandomSeed != 0 )
       {
-      registrationMethod->MetricSamplingReinitializeSeed(
-        atoi( antsRandomSeed ) );
+      registrationMethod->MetricSamplingReinitializeSeed( this->m_RegistrationRandomSeed );
       }
 
     for( unsigned int n = 0; n < stageMetricList.size(); n++ )
@@ -1014,6 +1018,8 @@ class RegistrationHelper : public itk::Object
   RealType                                m_UpperQuantile;
   std::ostream *                          m_LogStream;
 
+  int m_RegistrationRandomSeed;
+
   bool         m_ApplyLinearTransformsToFixedImageHeader;
   unsigned int m_PrintSimilarityMeasureInterval;
   unsigned int m_WriteIntervalVolumes;

diff --git a/Examples/itkantsRegistrationHelper.hxx b/Examples/itkantsRegistrationHelper.hxx
@@ -1825,6 +1825,11 @@ RegistrationHelper<TComputeType, VImageDimension>
         typename DisplacementFieldRegistrationType::Pointer displacementFieldRegistration =
           DisplacementFieldRegistrationType::New();
 
+	if ( this->m_RegistrationRandomSeed != 0 )
+	  {
+	  displacementFieldRegistration->MetricSamplingReinitializeSeed( this->m_RegistrationRandomSeed );
+	  }
+
         if( this->m_RestrictDeformationOptimizerWeights.size() > currentStageNumber )
           {
           if( this->m_RestrictDeformationOptimizerWeights[currentStageNumber].size() == VImageDimension )

diff --git a/Scripts/antsRegistrationSyNQuick.sh b/Scripts/antsRegistrationSyNQuick.sh
@@ -175,6 +175,8 @@ Optional arguments:
 
      -z:  collapse output transforms (default = 1)
 
+     -e:  Fix random seed to an int value (default = system time)
+
      NB:  Multiple image pairs can be specified for registration during the SyN stage.
           Specify additional images using the '-m' and '-f' options.  Note that image
           pair correspondence is given by the order specified on the command line.
@@ -282,9 +284,10 @@ NUMBEROFBINS=32
 MASKIMAGES=()
 USEHISTOGRAMMATCHING=0
 COLLAPSEOUTPUTTRANSFORMS=1
+RANDOMSEED=0
 
 # reading command line arguments
-while getopts "d:f:h:i:m:j:n:o:p:r:s:t:x:z:" OPT
+while getopts "d:e:f:h:i:m:j:n:o:p:r:s:t:x:z:" OPT
   do
   case $OPT in
       h) #help
@@ -293,6 +296,9 @@ while getopts "d:f:h:i:m:j:n:o:p:r:s:t:x:z:" OPT
    ;;
       d)  # dimensions
    DIM=$OPTARG
+   ;;
+      e)  # seed
+   RANDOMSEED=$OPTARG
    ;;
       x)  # inclusive mask
    MASKIMAGES[${#MASKIMAGES[@]}]=$OPTARG
@@ -563,7 +569,13 @@ case "$PRECISIONTYPE" in
   ;;
 esac
 
-COMMAND="${ANTS} --verbose 1 \
+RANDOMOPT=""
+
+if [[ ! $RANDOMSEED -eq 0 ]]; then
+    RANDOMOPT=" --random-seed $RANDOMSEED "
+fi
+
+COMMAND="${ANTS} --verbose 1 $RANDOMOPT \
                  --dimensionality $DIM $PRECISION \
                  --collapse-output-transforms $COLLAPSEOUTPUTTRANSFORMS \
                  --output [$OUTPUTNAME,${OUTPUTNAME}Warped.nii.gz,${OUTPUTNAME}InverseWarped.nii.gz] \