Vector Migration by Genes. Now rebased off of latest.

Eradication/GeneticProbability.cpp

@@ -242,6 +242,32 @@ namespace Kernel
         return !operator==( rThat );
     }
 
+    bool AlleleCombo::Compare( const AlleleCombo& rLeft, const AlleleCombo& rRight )
+    {
+        if( rLeft.GetNumLoci() < rRight.GetNumLoci() )
+        {
+            return true;
+        }
+        else if( rLeft.GetNumLoci() > rRight.GetNumLoci() )
+        {
+            return false;
+        }
+        else
+        {
+            if( rLeft.GetNumPossibleGenomes() == rRight.GetNumPossibleGenomes() )
+            {
+                // If they have the same number of possible genomes, we need some way to sort them so that the order
+                // is more consistent.  For example, we always want (a1,*) to come after (a0,*).
+                return ( rLeft.GetPossibleGenome( 0 ) < rRight.GetPossibleGenome( 0 ) );
+            }
+            else
+            {
+                // more possible genomes means more options so less specificity
+                return ( rLeft.GetNumPossibleGenomes() > rRight.GetNumPossibleGenomes() );
+            }
+        }
+    }
+
     bool AlleleCombo::IsNull() const
     {
         return (m_pAlleleCombo == nullptr);
@@ -1216,28 +1242,7 @@ namespace Kernel
 
     bool GeneticProbability::compareACP( const AlleleComboProbability& rLeft, const AlleleComboProbability& rRight )
     {
-        if( rLeft.GetAlleleCombo().GetNumLoci() < rRight.GetAlleleCombo().GetNumLoci() )
-        {
-            return true;
-        }
-        else if( rLeft.GetAlleleCombo().GetNumLoci() > rRight.GetAlleleCombo().GetNumLoci() )
-        {
-            return false;
-        }
-        else
-        {
-            if( rLeft.GetAlleleCombo().GetNumPossibleGenomes() == rRight.GetAlleleCombo().GetNumPossibleGenomes() )
-            {
-                // If they have the same number of possible genomes, we need some way to sort them so that the orde
-                // is more consistent.  For example, we always want (a1,*) to come after (a0,*).
-                return (rLeft.GetAlleleCombo().GetPossibleGenome( 0 ) < rRight.GetAlleleCombo().GetPossibleGenome( 0 ));
-            }
-            else
-            {
-                // more possible genomes means more options so less specificity
-                return (rLeft.GetAlleleCombo().GetNumPossibleGenomes() > rRight.GetAlleleCombo().GetNumPossibleGenomes());
-            }
-        }
+        return rLeft.GetAlleleCombo().Compare( rLeft.GetAlleleCombo(), rRight.GetAlleleCombo());
     }
 
     float GeneticProbability::AdditionOperation( float lhs, float rhs )

Eradication/GeneticProbability.h

@@ -63,6 +63,7 @@ namespace Kernel
 
         bool operator==( const AlleleCombo& rThat ) const;
         bool operator!=( const AlleleCombo& rThat ) const;
+        static bool Compare( const AlleleCombo& rLeft, const AlleleCombo& rRight );
 
         bool IsNull() const;
         int GetSpeciesIndex() const;
@@ -84,7 +85,7 @@ namespace Kernel
     };
 
     // AlleleComboProbability is the combination of an AlleleCombo and a single probability value.
-    // This is basically the collction of genomes with this specific probability.
+    // This is basically the collection of genomes with this specific probability.
     class AlleleComboProbability
     {
     public:

Eradication/IMigrationInfoVector.h

@@ -5,6 +5,7 @@
 #include "IMigrationInfo.h"
 #include "SimulationEnums.h"
 #include "VectorEnums.h"
+#include "VectorSpeciesParameters.h"
 
 
 namespace Kernel
@@ -22,6 +23,10 @@ namespace Kernel
                                                        const std::string& rSpeciesID,
                                                        IVectorSimulationContext* pivsc ) = 0;
         virtual Gender::Enum              ConvertVectorGender( VectorGender::Enum vector_gender ) const = 0;
+        virtual const std::vector<float>& GetFractionTraveling( VectorGender::Enum vector_gender, int by_genome_index ) = 0;
+        virtual int                       GetMigrationDataIndex( int species_index, VectorGenome& rGenome ) = 0;
+        virtual bool                      CanTravel() = 0;
+        virtual bool                      TravelByAlleles() = 0;
     };
 
     struct IMigrationInfoFactoryVector
@@ -32,6 +37,7 @@ namespace Kernel
                                                          const ::Configuration* config ) = 0;
         virtual IMigrationInfoVector* CreateMigrationInfoVector( const std::string& idreference,
                                                                  INodeContext *parent_node, 
-                                                                 const boost::bimap<ExternalNodeId_t, suids::suid>& rNodeIdSuidMap ) = 0;
+                                                                 const boost::bimap<ExternalNodeId_t, suids::suid>& rNodeIdSuidMap,
+                                                                 const VectorSpeciesParameters* m_species_params ) = 0;
     };
 }

Eradication/Insecticides.cpp

@@ -74,7 +74,7 @@ namespace Kernel
         float modifier_blocking = 1.0f;
         float modifier_killing = 1.0f;
 
-        initConfigTypeMap( "Allele_Combinations", &combo_strings, Insecticide_Allele_Combinations_DESC_TEXT, constraint_schema, allowed_values );
+        initConfigTypeMap( "Allele_Combinations",     &combo_strings,      Insecticide_Allele_Combinations_DESC_TEXT,     constraint_schema, allowed_values );
         initConfigTypeMap( "Larval_Killing_Modifier", &modifier_larval,    Insecticide_Larval_Killing_Modifier_DESC_TEXT, 0.0f, FLT_MAX, 1.0f );
         initConfigTypeMap( "Repelling_Modifier",      &modifier_repelling, Insecticide_Repelling_Modifier_DESC_TEXT,      0.0f, FLT_MAX, 1.0f );
         initConfigTypeMap( "Blocking_Modifier",       &modifier_blocking,  Insecticide_Blocking_Modifier_DESC_TEXT,       0.0f, FLT_MAX, 1.0f );

Eradication/Migration.cpp

@@ -75,7 +75,7 @@ namespace Kernel
         // --- Ages that are outside the range given in the metadata file are capped
         // --- at either the first or last value regardless if it is interpolation or piecewise
         // ------------------------------------------------------------------------------------
-        float rate; // = 0.0;
+        float rate = 0; 
         if( ageYears <= first_age_years )
         {
             rate = first_rate;
@@ -179,6 +179,7 @@ namespace Kernel
     void MigrationInfoFixedRate::Initialize( const std::vector<std::vector<MigrationRateData>>& rRateData )
     {
         m_RateCDF.clear();
+        m_TotalRate = 0.0;
         if( rRateData.size() > 0 )
         {
             for( auto& mrd : rRateData[0] )
@@ -189,6 +190,7 @@ namespace Kernel
                 release_assert( mrd.GetNumRates() == 1 );
 
                 m_RateCDF.push_back( mrd.GetRate( 0.0 ) );
+                m_TotalRate += mrd.GetRate( 0.0 );
             }
             NormalizeRates( m_RateCDF, m_TotalRate );
         }
@@ -265,27 +267,26 @@ namespace Kernel
 
     void MigrationInfoFixedRate::NormalizeRates( std::vector<float>& r_rate_cdf, float& r_total_rate )
     {
-        //  Calculate total migration rate
-        r_total_rate = 0.0;
-        for( int i = 0; i < r_rate_cdf.size(); i++)
-        {
-            r_total_rate += r_rate_cdf[i];
-        }
-
         if( (r_rate_cdf.size() > 0) && (r_total_rate > 0.0) )
         {
             //  Set probability of each location
+            // Values in the migration_rate_cdf[] are compared against values generated by RANDOMBASE::e(), which 
+            // is guaranteed to create values between 0.0f and 1.0f.  We need to explicitly set values that are 
+            // floating point rounding error close to 1 to 1, otherwise floating-point rounding errors might result 
+            // in running past the end of valid values in the array when picking a migration destination
             r_rate_cdf[0] /= r_total_rate;
             for(int i = 1; i < r_rate_cdf.size(); i++)
             {
-                r_rate_cdf[i] = (r_rate_cdf[i] / r_total_rate) + r_rate_cdf[i-1];
+                float tmp_rate = ( r_rate_cdf[i] / r_total_rate ) + r_rate_cdf[i - 1];
+                if( 1 - tmp_rate < FLT_EPSILON ) // if close to 1, call it 1
+                {
+                    r_rate_cdf[i] = 1;
+                }
+                else
+                {
+                    r_rate_cdf[i] = tmp_rate;
+                }
             }
-
-            // Values in the migration_rate_cdf[] are compared against values generated by RANDOMBASE::e(), which 
-            // is guaranteed to create values between 0.0f and 1.0f.  We need to explicitly set the last value in
-            // migration_rate_cdf to 1.0f, otherwise floating-point rounding errors might result in running past the
-            // end of valid values in the array when picking a migration destination
-            r_rate_cdf[ r_rate_cdf.size() - 1 ] = 1.0f;
         }
     }
 
@@ -355,6 +356,7 @@ namespace Kernel
         {
             float rate = mrd.GetRate( ageYears );
             m_RateCDF.push_back( rate );
+            m_TotalRate += rate;
         }
 
         SaveRawRates(m_RateCDF, gender);
@@ -541,7 +543,7 @@ namespace Kernel
                                 rRateData[ ig ].push_back( mrd );
                             }
                             else if( (rRateData[ ig ][ node_index ].GetToNodeSuid()      != to_node_id_suid) ||
-                                     (rRateData[ ig ][ node_index ].GetMigrationType() != m_MigrationType) )
+                                     (rRateData[ ig ][ node_index ].GetMigrationType()   != m_MigrationType) )
                             {
                                 std::stringstream ss;
                                 ss << "In file '" << m_Filename << "', the 'To' Node IDs are not the same for the Age Data sections for fromNodeId = " << fromNodeID;
@@ -642,7 +644,8 @@ static const char* NODE_OFFSETS          = "NodeOffsets";            // required
                     std::stringstream ss;
                     ss << metadata_filepath << "[" << METADATA << "][" << MD_MIGRATION_TYPE << "]";
                     throw IncoherentConfigurationException( __FILE__, __LINE__, __FUNCTION__, 
-                                                            "m_MigrationType", exp_mig_type_str.c_str(),
+                                                            "m_MigrationType", 
+                                                            exp_mig_type_str.c_str(),
                                                             ss.str().c_str(),
                                                             file_mig_type_str.c_str() );
                 }
@@ -660,6 +663,16 @@ static const char* NODE_OFFSETS          = "NodeOffsets";            // required
                 {
                     THROW_ENUM_EXCEPTION( GenderDataType, metadata_filepath, MD_GENDER_DATA_TYPE, gd_str )
                 }
+                else if( m_GenderDataType == GenderDataType::VECTOR_MIGRATION_BY_GENETICS )
+                {
+                    std::stringstream ss;
+                    ss << metadata_filepath << "[" << METADATA << "][" << MD_GENDER_DATA_TYPE << "] ";
+                    throw IncoherentConfigurationException( __FILE__, __LINE__, __FUNCTION__,
+                                                            "GenderDataType",
+                                                            gd_str.c_str(),
+                                                            ss.str().c_str(),
+                                                            "cannot use for human migration");
+                }
             }
 
             // ---------------------------------------------------------------------------------
@@ -810,9 +823,9 @@ static const char* NODE_OFFSETS          = "NodeOffsets";            // required
     uint32_t MigrationInfoFile::GetNumGenderDataChunks() const
     {
         uint32_t num_gender_data_chunks; // = 1;
-        if( m_GenderDataType == GenderDataType::SAME_FOR_BOTH_GENDERS )
+        if( m_GenderDataType != GenderDataType::ONE_FOR_EACH_GENDER )
         {
-            num_gender_data_chunks = 1; // genders use same data
+            num_gender_data_chunks = 1; // genders use same data, including VECTOR_MIGRATION_BY_GENETICS
         }
         else
         {

Eradication/Migration.h

@@ -28,8 +28,10 @@ namespace std
 namespace Kernel
 {
     ENUM_DEFINE(GenderDataType,
-        ENUM_VALUE_SPEC(SAME_FOR_BOTH_GENDERS , 0)  // The one set of the data is used for both genders
-        ENUM_VALUE_SPEC(ONE_FOR_EACH_GENDER   , 1)) // There are two sets of data - one for each gender
+        ENUM_VALUE_SPEC( SAME_FOR_BOTH_GENDERS       , 0)   // The one set of the data is used for both genders
+        ENUM_VALUE_SPEC( ONE_FOR_EACH_GENDER         , 1)   // There are two sets of data - one for each gender
+        ENUM_VALUE_SPEC( VECTOR_MIGRATION_BY_GENETICS, 2 )) // Only for vector migration: we will be using 
+                                                            // vector genetics (which accounts for gender) for migration
 
     ENUM_DEFINE(InterpolationType,
         ENUM_VALUE_SPEC(LINEAR_INTERPOLATION , 0)  // Interpolate between ages - no extrapolation
@@ -84,10 +86,10 @@ namespace Kernel
                                         MigrationType::Enum &migration_type,
                                         float &time, 
                                         float dt = FLT_MAX ) override; // FLT_MAX for humans, dt for vectors
-        virtual void SetContextTo( INodeContext* _parent ) override;
+        virtual void  SetContextTo( INodeContext* _parent ) override;
         virtual float GetTotalRate( Gender::Enum gender = Gender::MALE ) const override;
-        virtual const std::vector<float>& GetCumulativeDistributionFunction( Gender::Enum gender = Gender::MALE ) const override;
-        virtual const std::vector<suids::suid>& GetReachableNodes( Gender::Enum gender = Gender::MALE ) const override;
+        virtual const std::vector<float>&               GetCumulativeDistributionFunction( Gender::Enum gender = Gender::MALE ) const override;
+        virtual const std::vector<suids::suid>&         GetReachableNodes( Gender::Enum gender = Gender::MALE ) const override;
         virtual const std::vector<MigrationType::Enum>& GetMigrationTypes( Gender::Enum gender = Gender::MALE ) const override;
 
         virtual bool IsHeterogeneityEnabled() const override;
@@ -215,7 +217,9 @@ namespace Kernel
                                const boost::bimap<ExternalNodeId_t, suids::suid>& rNodeidSuidMap,
                                std::vector<std::vector<MigrationRateData>>&       rRateData );
 
-        MigrationType::Enum GetMigrationType() const { return m_MigrationType; }
+        MigrationType::Enum  GetMigrationType()  const { return m_MigrationType; }
+        GenderDataType::Enum GetGenderDataType() const { return m_GenderDataType; }
+
 
     protected:
         // Returns the expected size of the binary file