diff --git a/CHANGES.html b/CHANGES.html
index 3928c8df4..1b800e4c9 100644
--- a/CHANGES.html
+++ b/CHANGES.html
@@ -56,11 +56,15 @@ us a note on ns-developers mailing list.</p>
 <ul>
 <li>The helpers of the NetDevices supporting flow control (PointToPointHelper, CsmaHelper, SimpleNetDeviceHelper, WifiHelper) now provide a <b>DisableFlowControl</b> method to disable flow control. If flow control is disabled, the Traffic Control layer forwards packets down to the NetDevice even if there is no room for them in the NetDevice queue(s)</li>
 <li>Added a new trace source <b>TcDrop</b> in TrafficControlLayer for tracing packets that have been dropped because no queue disc is installed on the device, the device supports flow control and the device queue is full.</li>
+<li>Added a new class <b>PhasedArraySpectrumPropagationLossModel</b>, and its <b>DoCalcRxPowerSpectralDensity</b> function has two additional parameters: TX and RX antenna arrays. Should be inherited by models that need to know antenna arrays in order to calculate RX PSD.</li>
 </ul>
 <h2>Changes to existing API:</h2>
 <ul>
 <li>internet: Support for Network Simulation Cradle (NSC) TCP has been removed.</li>
 <li>fd-net-device: Support for PlanetLabFdNetDeviceHelper has been removed.</li>
+<li><b>ThreeGppSpectrumPropagationLossModel</b> now inherits <b>PhasedArraySpectrumPropagationLossModel</b>. The modules that use <b>ThreeGppSpectrumPropagationLossModel</b> should implement <b>SpectrumPhy::GetAntenna</b> that will return the instance of <b>PhasedArrayModel</b>.</li>
+<li><b>AddDevice</b> function is removed from <b>ThreeGppSpectrumPropagationLossModel</b> to support multiple arrays per device.</li>
+<li><b>SpectrumPhy</b> function <b>GetRxAntenna</b> is renamed to <b>GetAntenna</b>, and its return value is changed to <b>Ptr<Object></b> instead of <b>Ptr<AntennaModel></b> to support also <b>PhasedArrayModel</b> type of antenna.</li>
 </ul>
 <h2>Changes to build system:</h2>
 <ul>
diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md
index 0c49654e1..0426f9e1b 100644
--- a/RELEASE_NOTES.md
+++ b/RELEASE_NOTES.md
@@ -26,8 +26,11 @@ requirements (Note: not all ns-3 features are available on all systems):
 
 ### New user-visible features
 
+- (spectrum) ThreeGppSpectrumPropagationLossModel and ThreeGppChannelModel now support multiple PhasedArrayModel instances per device. This feature can be used to implement MIMO.
 ### Bugs fixed
 
+- (spectrum) Fix condition for channel matrix update in ThreeGppChannelModel (left and right operand were pointing to the same object)
+- (spectrum) Assign stream to random variable for doppler term in ThreeGppChannelModel (moved from ThreeGppSpectrumPropagationLossModel)
 - (wifi) #467 - WiFi: Failed association process
 - (wifi) #468 - WiFi: Wrong txDuration for trigger frame
 - (wifi) #475 - Wi-Fi: Assert when sending OFDMA DL to STAs with different TIDs
diff --git a/examples/channel-models/three-gpp-v2v-channel-example.cc b/examples/channel-models/three-gpp-v2v-channel-example.cc
index cb0b10c9b..15ab147eb 100644
--- a/examples/channel-models/three-gpp-v2v-channel-example.cc
+++ b/examples/channel-models/three-gpp-v2v-channel-example.cc
@@ -54,6 +54,42 @@ static Ptr<ThreeGppPropagationLossModel> m_propagationLossModel; //!< the Propag
 static Ptr<ThreeGppSpectrumPropagationLossModel> m_spectrumLossModel; //!< the SpectrumPropagationLossModel object
 static Ptr<ChannelConditionModel> m_condModel; //!< the ChannelConditionModel object
 
+/*
+ * \brief A structure that holds the parameters for the ComputeSnr
+ * function. In this way the problem with the limited
+ * number of parameters of method Schedule is avoided.
+ */
+struct ComputeSnrParams
+{
+  Ptr<MobilityModel> txMob; //!< the tx mobility model
+  Ptr<MobilityModel> rxMob; //!< the rx mobility model
+  Ptr<SpectrumValue> txPsd; //!< the PSD of the tx signal
+  double noiseFigure; //!< the noise figure in dB
+  Ptr<PhasedArrayModel> txAntenna; //!< the tx antenna array
+  Ptr<PhasedArrayModel> rxAntenna; //!< the rx antenna array
+
+  /**
+   * \brief Constructor
+   * \param pTxMob the tx mobility model
+   * \param pRxMob the rx mobility model
+   * \param pTxPsd the PSD of the tx signal
+   * \param pNoiseFigure the noise figure in dB
+   * \param pTxAntenna the tx antenna array
+   * \param pRxAntenna the rx antenna array
+   */
+  ComputeSnrParams (Ptr<MobilityModel> pTxMob, Ptr<MobilityModel> pRxMob,
+                    Ptr<SpectrumValue> pTxPsd, double pNoiseFigure,
+                    Ptr<PhasedArrayModel> pTxAntenna, Ptr<PhasedArrayModel> pRxAntenna)
+  {
+    txMob = pTxMob;
+    rxMob = pRxMob;
+    txPsd = pTxPsd;
+    noiseFigure = pNoiseFigure;
+    txAntenna = pTxAntenna;
+    rxAntenna = pRxAntenna;
+  }
+};
+
 /**
  * Perform the beamforming using the DFT beamforming method
  * \param thisDevice the device performing the beamforming
@@ -78,36 +114,33 @@ DoBeamforming (Ptr<NetDevice> thisDevice, Ptr<PhasedArrayModel> thisAntenna, Ptr
 
 /**
  * Compute the average SNR
- * \param txMob the tx mobility model
- * \param rxMob the rx mobility model
- * \param txPsd the PSD of the transmitting signal
- * \param noiseFigure the noise figure in dB
+ * \param params A structure that holds a bunch of parameters needed by ComputSnr function to calculate the average SNR
  */
 static void
-ComputeSnr (Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, Ptr<const SpectrumValue> txPsd, double noiseFigure)
+ComputeSnr (ComputeSnrParams& params)
 {
-  Ptr<SpectrumValue> rxPsd = txPsd->Copy ();
+  Ptr<SpectrumValue> rxPsd = params.txPsd->Copy ();
 
   // check the channel condition
-  Ptr<ChannelCondition> cond = m_condModel->GetChannelCondition (txMob, rxMob);
+  Ptr<ChannelCondition> cond = m_condModel->GetChannelCondition (params.txMob, params.rxMob);
 
   // apply the pathloss
-  double propagationGainDb = m_propagationLossModel->CalcRxPower (0, txMob, rxMob);
+  double propagationGainDb = m_propagationLossModel->CalcRxPower (0, params.txMob, params.rxMob);
   NS_LOG_DEBUG ("Pathloss " << -propagationGainDb << " dB");
   double propagationGainLinear = std::pow (10.0, (propagationGainDb) / 10.0);
   *(rxPsd) *= propagationGainLinear;
 
   // apply the fast fading and the beamforming gain
-  rxPsd = m_spectrumLossModel->CalcRxPowerSpectralDensity (rxPsd, txMob, rxMob);
+  rxPsd = m_spectrumLossModel->CalcRxPowerSpectralDensity (rxPsd, params.txMob, params.rxMob, params.txAntenna, params.rxAntenna);
   NS_LOG_DEBUG ("Average rx power " << 10 * log10 (Sum (*rxPsd) * 180e3) << " dB");
 
   // create the noise psd
   // taken from lte-spectrum-value-helper
   const double kT_dBm_Hz = -174.0; // dBm/Hz
   double kT_W_Hz = std::pow (10.0, (kT_dBm_Hz - 30) / 10.0);
-  double noiseFigureLinear = std::pow (10.0, noiseFigure / 10.0);
+  double noiseFigureLinear = std::pow (10.0, params.noiseFigure / 10.0);
   double noisePowerSpectralDensity =  kT_W_Hz * noiseFigureLinear;
-  Ptr<SpectrumValue> noisePsd = Create <SpectrumValue> (txPsd->GetSpectrumModel ());
+  Ptr<SpectrumValue> noisePsd = Create <SpectrumValue> (params.txPsd->GetSpectrumModel ());
   (*noisePsd) = noisePowerSpectralDensity;
 
   // compute the SNR
@@ -117,10 +150,10 @@ ComputeSnr (Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, Ptr<const Spectr
   std::ofstream f;
   f.open ("example-output.txt", std::ios::out | std::ios::app);
   f << Simulator::Now ().GetSeconds () << " " // time [s]
-    << txMob->GetPosition ().x << " "
-    << txMob->GetPosition ().y << " "
-    << rxMob->GetPosition ().x << " "
-    << rxMob->GetPosition ().y << " "
+    << params.txMob->GetPosition ().x << " "
+    << params.txMob->GetPosition ().y << " "
+    << params.rxMob->GetPosition ().x << " "
+    << params.rxMob->GetPosition ().y << " "
     << cond->GetLosCondition () << " " // channel state
     << 10 * log10 (Sum (*rxPsd) / Sum (*noisePsd)) << " " // SNR [dB]
     << -propagationGainDb << std::endl; // pathloss [dB]
@@ -301,10 +334,6 @@ main (int argc, char *argv[])
   m_spectrumLossModel = CreateObjectWithAttributes<ThreeGppSpectrumPropagationLossModel> ("ChannelModel", PointerValue (channelModel));
   m_spectrumLossModel->SetAttribute ("vScatt", DoubleValue (vScatt));
 
-  // initialize the devices in the ThreeGppSpectrumPropagationLossModel
-  m_spectrumLossModel->AddDevice (txDev, txAntenna);
-  m_spectrumLossModel->AddDevice (rxDev, rxAntenna);
-
   BuildingsHelper::Install (nodes);
 
   // set the beamforming vectors
@@ -332,7 +361,7 @@ main (int argc, char *argv[])
 
   for (int i = 0; i < simTime / timeRes; i++)
     {
-      Simulator::Schedule (timeRes * i, &ComputeSnr, txMob, rxMob, txPsd, noiseFigure);
+      Simulator::Schedule (timeRes * i, &ComputeSnr, ComputeSnrParams (txMob, rxMob, txPsd, noiseFigure, txAntenna, rxAntenna));
     }
 
   // initialize the output file
diff --git a/src/spectrum/CMakeLists.txt b/src/spectrum/CMakeLists.txt
index 072952b12..72e387b2a 100644
--- a/src/spectrum/CMakeLists.txt
+++ b/src/spectrum/CMakeLists.txt
@@ -27,6 +27,7 @@ set(source_files
     model/spectrum-model.cc
     model/spectrum-phy.cc
     model/spectrum-propagation-loss-model.cc
+    model/phased-array-spectrum-propagation-loss-model.cc
     model/spectrum-signal-parameters.cc
     model/spectrum-value.cc
     model/three-gpp-channel-model.cc
@@ -64,6 +65,7 @@ set(header_files
     model/spectrum-model.h
     model/spectrum-phy.h
     model/spectrum-propagation-loss-model.h
+    model/phased-array-spectrum-propagation-loss-model.h
     model/spectrum-signal-parameters.h
     model/spectrum-value.h
     model/three-gpp-channel-model.h
diff --git a/src/spectrum/doc/spectrum.rst b/src/spectrum/doc/spectrum.rst
index cba0ce038..b97912665 100644
--- a/src/spectrum/doc/spectrum.rst
+++ b/src/spectrum/doc/spectrum.rst
@@ -123,7 +123,7 @@ The module provides two ``SpectrumChannel`` implementations:
 ``SingleModelSpectrumChannel`` and ``MultiModelSpectrumChannel``. They
 both provide this functionality:
 
- * Propagation loss modeling, in two forms:
+ * Propagation loss modeling, in three forms:
 
    - you can plug models based on ``PropagationLossModel`` on these
      channels. Only linear models (where the loss value does not
@@ -136,6 +136,14 @@ both provide this functionality:
      a separate loss value is calculated and applied to each component
      of the power spectral density.
 
+   - you can plug models based on ``PhasedArraySpectrumPropagationLossModel`` 
+     on these channels. These models can have frequency-dependent loss, i.e.,
+     a separate loss value is calculated and applied to each component
+     of the power spectral density. Additionally, these models support 
+     the phased antenna array at the transmitter and the receiver, i.e., 
+     ns-3 antenna type ``PhasedArrayModel``.
+
+
  * Propagation delay modeling, by plugging a model based on
    ``PropagationDelayModel``. The delay is independent of frequency and
    applied to the signal as a whole. Delay modeling is implemented by
@@ -642,26 +650,24 @@ ThreeGppChannelModel and ThreeGppSpectrumPropagationLossModel.
 ThreeGppSpectrumPropagationLossModel
 ####################################
 
-The class ThreeGppSpectrumPropagationLossModel extends the SpectrumPropagationLossModel
-interface and enables the modeling of frequency
-dependent propagation phenomena. The main method is DoCalcRxPowerSpectralDensity,
-which takes as input the power spectral density (PSD) of the transmitted signal,
-the mobility models of the transmitting node and receiving node, and
-returns the PSD of the received signal.
+The class ThreeGppSpectrumPropagationLossModel implements the 
+PhasedArraySpectrumPropagationLossModel interface and enables the modeling of frequency
+dependent propagation phenomena while taking into account the specific pair of the 
+antenna array at the transmitter and the receiver. The main method is 
+DoCalcRxPowerSpectralDensity, which takes as input the power spectral density (PSD) 
+of the transmitted signal, the mobility models of the transmitting node and receiving node, 
+and the antenna array of the transmitting node, and of the receiving node. 
+Finally, it returns the PSD of the received signal.
 
 Procedure used to compute the PSD of to compute the PSD of the received signal:
 
 1. Retrieve the beamforming vectors
 To account for the beamforming, ThreeGppSpectrumPropagationLossModel has to
 retrieve the beamforming vectors of the transmitting and receiving antennas.
-The method DoCalcRxPowerSpectralDensity uses m_deviceAntennaMap to obtain the
-antenna objects associated to the transmitting and receiving devices, and calls
-the method GetCurrentBeamformingVector to retrieve the beamforming vectors.
-For each device using the channel, the m_deviceAntennaMap contains the associated
-antenna object of type PhasedArrayModel. Since the mapping is one-to-one,
-the model supports a single antenna object for each device.
-The m_deviceAntennaMap has to be initialized by inserting the device-antenna
-pairs using the method AddDevice.
+The method DoCalcRxPowerSpectralDensity uses the antenna objects 
+that are passed as parameters for both the transmitting and receiving devices, 
+and calls the method GetCurrentBeamformingVector to retrieve the beamforming vectors 
+of these antenna objects.
 
 2. Retrieve the channel matrix
 The ThreeGppSpectrumPropagationLossModel relies on the ThreeGppChannelModel class
diff --git a/src/spectrum/examples/three-gpp-channel-example.cc b/src/spectrum/examples/three-gpp-channel-example.cc
index bbeb6c2eb..12cb5bac8 100644
--- a/src/spectrum/examples/three-gpp-channel-example.cc
+++ b/src/spectrum/examples/three-gpp-channel-example.cc
@@ -50,6 +50,41 @@ using namespace ns3;
 static Ptr<ThreeGppPropagationLossModel> m_propagationLossModel; //!< the PropagationLossModel object
 static Ptr<ThreeGppSpectrumPropagationLossModel> m_spectrumLossModel; //!< the SpectrumPropagationLossModel object
 
+/**
+ * \brief A structure that holds the parameters for the
+ * ComputeSnr function. In this way the problem with the limited
+ * number of parameters of method Schedule is avoided.
+ */
+struct ComputeSnrParams
+{
+  Ptr<MobilityModel> txMob; //!< the tx mobility model
+  Ptr<MobilityModel> rxMob; //!< the rx mobility model
+  double txPow; //!< the tx power in dBm
+  double noiseFigure; //!< the noise figure in dB
+  Ptr<PhasedArrayModel> txAntenna; //!< the tx antenna array
+  Ptr<PhasedArrayModel> rxAntenna; //!< the rx antenna array
+
+  /**
+   * \brief Constructor
+   * \param pTxMob the tx mobility model
+   * \param pRxMob the rx mobility model
+   * \param pTxPow the tx power in dBm
+   * \param pNoiseFigure the noise figure in dB
+   * \param pTxAntenna the tx antenna array
+   * \param pRxAntenna the rx antenna array
+   */
+  ComputeSnrParams (Ptr<MobilityModel> pTxMob, Ptr<MobilityModel> pRxMob, double pTxPow, double pNoiseFigure,
+                    Ptr<PhasedArrayModel> pTxAntenna, Ptr<PhasedArrayModel> pRxAntenna)
+  {
+    txMob = pTxMob;
+    rxMob = pRxMob;
+    txPow = pTxPow;
+    noiseFigure = pNoiseFigure;
+    txAntenna = pTxAntenna;
+    rxAntenna = pRxAntenna;
+  }
+};
+
 /**
  * Perform the beamforming using the DFT beamforming method
  * \param thisDevice the device performing the beamforming
@@ -93,13 +128,10 @@ DoBeamforming (Ptr<NetDevice> thisDevice, Ptr<PhasedArrayModel> thisAntenna, Ptr
 
 /**
  * Compute the average SNR
- * \param txMob the tx mobility model
- * \param rxMob the rx mobility model
- * \param txPow the transmitting power in dBm
- * \param noiseFigure the noise figure in dB
+ * \param params A structure that holds the parameters that are needed to perform calculations in ComputeSnr
  */
 static void
-ComputeSnr (Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, double txPow, double noiseFigure)
+ComputeSnr (ComputeSnrParams& params)
 {
   // Create the tx PSD using the LteSpectrumValueHelper
   // 100 RBs corresponds to 18 MHz (1 RB = 180 kHz)
@@ -109,23 +141,39 @@ ComputeSnr (Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, double txPow, do
   {
     activeRbs0[i] = i;
   }
-  Ptr<SpectrumValue> txPsd = LteSpectrumValueHelper::CreateTxPowerSpectralDensity (2100, 100, txPow, activeRbs0);
+  Ptr<SpectrumValue> txPsd = LteSpectrumValueHelper::CreateTxPowerSpectralDensity (2100, 100, params.txPow, activeRbs0);
   Ptr<SpectrumValue> rxPsd = txPsd->Copy ();
   NS_LOG_DEBUG ("Average tx power " << 10*log10(Sum (*txPsd) * 180e3) << " dB");
 
   // create the noise PSD
-  Ptr<SpectrumValue> noisePsd = LteSpectrumValueHelper::CreateNoisePowerSpectralDensity (2100, 100, noiseFigure);
+  Ptr<SpectrumValue> noisePsd = LteSpectrumValueHelper::CreateNoisePowerSpectralDensity (2100, 100, params.noiseFigure);
   NS_LOG_DEBUG ("Average noise power " << 10*log10 (Sum (*noisePsd) * 180e3) << " dB");
 
   // apply the pathloss
-  double propagationGainDb = m_propagationLossModel->CalcRxPower (0, txMob, rxMob);
+  double propagationGainDb = m_propagationLossModel->CalcRxPower (0, params.txMob, params.rxMob);
   NS_LOG_DEBUG ("Pathloss " << -propagationGainDb << " dB");
   double propagationGainLinear = std::pow (10.0, (propagationGainDb) / 10.0);
   *(rxPsd) *= propagationGainLinear;
 
+  if (params.txAntenna == nullptr)
+    {
+      NS_ABORT_MSG("params.txAntenna is nullptr!");
+    }
+  else if (params.rxAntenna == nullptr)
+    {
+      NS_ABORT_MSG("params.rxAntenna is nullptr!");
+    }
+
+
+  Ptr<const PhasedArrayModel> a = params.txAntenna;
+  if (a == nullptr)
+    {
+      NS_ABORT_MSG("params.txAntenna is nullptr!");
+    }
+
   // apply the fast fading and the beamforming gain
-  rxPsd = m_spectrumLossModel->CalcRxPowerSpectralDensity (rxPsd, txMob, rxMob);
-  NS_LOG_DEBUG ("Average rx power " << 10*log10 (Sum (*rxPsd) * 180e3) << " dB");
+  rxPsd = m_spectrumLossModel->CalcRxPowerSpectralDensity (rxPsd, params.txMob, params.rxMob, params.txAntenna, params.rxAntenna);
+  NS_LOG_DEBUG ("Average rx power " << 10 * log10 (Sum (*rxPsd) * 180e3) << " dB");
 
   // compute the SNR
   NS_LOG_DEBUG ("Average SNR " << 10 * log10 (Sum (*rxPsd) / Sum (*noisePsd)) << " dB");
@@ -231,17 +279,13 @@ main (int argc, char *argv[])
   Ptr<PhasedArrayModel> txAntenna = CreateObjectWithAttributes<UniformPlanarArray> ("NumColumns", UintegerValue (2), "NumRows", UintegerValue (2));
   Ptr<PhasedArrayModel> rxAntenna = CreateObjectWithAttributes<UniformPlanarArray> ("NumColumns", UintegerValue (2), "NumRows", UintegerValue (2));
 
-  // initialize the devices in the ThreeGppSpectrumPropagationLossModel
-  m_spectrumLossModel->AddDevice (txDev, txAntenna);
-  m_spectrumLossModel->AddDevice (rxDev, rxAntenna);
-
   // set the beamforming vectors
   DoBeamforming (txDev, txAntenna, rxDev);
   DoBeamforming (rxDev, rxAntenna, txDev);
 
   for (int i = 0; i < floor (simTime / timeRes); i++)
   {
-    Simulator::Schedule (MilliSeconds (timeRes*i), &ComputeSnr, txMob, rxMob, txPow, noiseFigure);
+    Simulator::Schedule (MilliSeconds (timeRes*i), &ComputeSnr, ComputeSnrParams (txMob, rxMob, txPow, noiseFigure, txAntenna, rxAntenna));
   }
 
   Simulator::Run ();
diff --git a/src/spectrum/model/matrix-based-channel-model.h b/src/spectrum/model/matrix-based-channel-model.h
index 37bccff96..85f6ff63a 100644
--- a/src/spectrum/model/matrix-based-channel-model.h
+++ b/src/spectrum/model/matrix-based-channel-model.h
@@ -112,14 +112,15 @@ public:
                                                Ptr<const PhasedArrayModel> bAntenna) = 0;
   
   /**
-   * Calculate the channel key using the Cantor function
-   * \param x1 first value
-   * \param x2 second value
-   * \return \f$ (((x1 + x2) * (x1 + x2 + 1))/2) + x2; \f$
+   * Generate a unique value for the pair of unsigned integer of 32 bits,
+   * where the order does not matter, i.e., the same value will be returned for (a,b) and (b,a).
+   * \param a the first value
+   * \param b the second value
+   * \return return an integer representing a unique value for the pair of values
    */
-  static constexpr uint32_t GetKey (uint32_t x1, uint32_t x2)
+  static uint64_t GetKey (uint32_t a, uint32_t b)
   {
-   return (((x1 + x2) * (x1 + x2 + 1)) / 2) + x2;
+   return (uint64_t) std::min (a, b) << 32 | std::max (a, b);
   }
 
   static const uint8_t AOA_INDEX = 0; //!< index of the AOA value in the m_angle array
diff --git a/src/spectrum/model/multi-model-spectrum-channel.cc b/src/spectrum/model/multi-model-spectrum-channel.cc
index de746186d..7017a26c6 100644
--- a/src/spectrum/model/multi-model-spectrum-channel.cc
+++ b/src/spectrum/model/multi-model-spectrum-channel.cc
@@ -269,6 +269,19 @@ MultiModelSpectrumChannel::StartTx (Ptr<SpectrumSignalParameters> txParams)
 
           if ((*rxPhyIterator) != txParams->txPhy)
             {
+              Ptr<NetDevice> rxNetDevice = (*rxPhyIterator)->GetDevice ();
+              Ptr<NetDevice> txNetDevice = txParams->txPhy->GetDevice ();
+
+              if (rxNetDevice && txNetDevice)
+                {
+                  // we assume that devices are attached to a node
+                  if (rxNetDevice->GetNode()->GetId() == txNetDevice->GetNode()->GetId())
+                    {
+                      NS_LOG_DEBUG ("Skipping the pathloss calculation among different antennas of the same node, not supported yet by any pathloss model in ns-3.");
+                      continue;
+                    }
+                }
+
               NS_LOG_LOGIC ("copying signal parameters " << txParams);
               Ptr<SpectrumSignalParameters> rxParams = txParams->Copy ();
               rxParams->psd = Copy<SpectrumValue> (convertedTxPowerSpectrum);
@@ -320,6 +333,21 @@ MultiModelSpectrumChannel::StartTx (Ptr<SpectrumSignalParameters> txParams)
                     {
                       rxParams->psd = m_spectrumPropagationLoss->CalcRxPowerSpectralDensity (rxParams->psd, txMobility, receiverMobility);
                     }
+                  else if (m_phasedArraySpectrumPropagationLoss)
+                    {
+                      // they cannot be combined, this case is not supported yet by any propagation model in ns-3
+                      NS_ASSERT_MSG (rxParams->txAntenna == nullptr && rxAntenna == nullptr, " Either AntennaModel or PhasedArrayModel can be used on a pair of TX/RX device");
+
+                      Ptr<PhasedArraySpectrumPhy> txPhasedArraySpectrumPhy = DynamicCast<PhasedArraySpectrumPhy> (txParams->txPhy);
+                      Ptr<PhasedArraySpectrumPhy> rxPhasedArraySpectrumPhy = DynamicCast<PhasedArraySpectrumPhy> (*rxPhyIterator);
+
+                      NS_ASSERT_MSG (txPhasedArraySpectrumPhy && rxPhasedArraySpectrumPhy, "PhasedArraySpectrumPhy should be installed at both TX and RX in order to use PhasedArraySpectrumPropagationLoss.");
+                      Ptr<const PhasedArrayModel> txPhasedArrayModel = txPhasedArraySpectrumPhy->GetPhasedArrayModel ();
+                      Ptr<const PhasedArrayModel> rxPhasedArrayModel = rxPhasedArraySpectrumPhy->GetPhasedArrayModel ();
+                      NS_ASSERT_MSG (txPhasedArrayModel && rxPhasedArrayModel, "PhasedArrayModel instances should be installed at both TX and RX SpectrumPhy in order to use PhasedArraySpectrumPropagationLoss.");
+
+                      rxParams->psd = m_phasedArraySpectrumPropagationLoss->CalcRxPowerSpectralDensity (rxParams->psd, txMobility, receiverMobility, txPhasedArrayModel, rxPhasedArrayModel);
+                     }
 
                   if (m_propagationDelay)
                     {
@@ -327,11 +355,10 @@ MultiModelSpectrumChannel::StartTx (Ptr<SpectrumSignalParameters> txParams)
                     }
                 }
 
-              Ptr<NetDevice> netDev = (*rxPhyIterator)->GetDevice ();
-              if (netDev)
+              if (rxNetDevice)
                 {
                   // the receiver has a NetDevice, so we expect that it is attached to a Node
-                  uint32_t dstNode =  netDev->GetNode ()->GetId ();
+                  uint32_t dstNode = rxNetDevice->GetNode ()->GetId ();
                   Simulator::ScheduleWithContext (dstNode, delay, &MultiModelSpectrumChannel::StartRx, this,
                                                   rxParams, *rxPhyIterator);
                 }
diff --git a/src/spectrum/model/phased-array-spectrum-propagation-loss-model.cc b/src/spectrum/model/phased-array-spectrum-propagation-loss-model.cc
new file mode 100644
index 000000000..703b1e2a1
--- /dev/null
+++ b/src/spectrum/model/phased-array-spectrum-propagation-loss-model.cc
@@ -0,0 +1,80 @@
+/* -*- Mode:C++; c-file-style:"gnu"; indent-tabs-mode:nil; -*- */
+/*
+ * Copyright (c) 2021 CTTC
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License version 2 as
+ * published by the Free Software Foundation;
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ *
+ */
+
+
+#include "phased-array-spectrum-propagation-loss-model.h"
+#include <ns3/log.h>
+#include <ns3/phased-array-model.h>
+
+namespace ns3 {
+
+NS_LOG_COMPONENT_DEFINE ("PhasedArraySpectrumPropagationLossModel");
+
+NS_OBJECT_ENSURE_REGISTERED (PhasedArraySpectrumPropagationLossModel);
+
+PhasedArraySpectrumPropagationLossModel::PhasedArraySpectrumPropagationLossModel ()
+  : m_next (0)
+{
+}
+
+PhasedArraySpectrumPropagationLossModel::~PhasedArraySpectrumPropagationLossModel ()
+{
+}
+
+void
+PhasedArraySpectrumPropagationLossModel::DoDispose ()
+{
+  m_next = 0;
+}
+
+TypeId
+PhasedArraySpectrumPropagationLossModel::GetTypeId (void)
+{
+  static TypeId tid = TypeId ("ns3::PhasedArraySpectrumPropagationLossModel")
+    .SetParent<Object> ()
+    .SetGroupName ("Spectrum")
+  ;
+  return tid;
+}
+
+
+void PhasedArraySpectrumPropagationLossModel::SetNext (Ptr<PhasedArraySpectrumPropagationLossModel> next)
+{
+  m_next = next;
+}
+
+Ptr<SpectrumValue>
+PhasedArraySpectrumPropagationLossModel::CalcRxPowerSpectralDensity (Ptr<const SpectrumValue> txPsd,
+                                                                     Ptr<const MobilityModel> a,
+                                                                     Ptr<const MobilityModel> b,
+                                                                     Ptr<const PhasedArrayModel> aPhasedArrayModel,
+                                                                     Ptr<const PhasedArrayModel> bPhasedArrayModel) const
+{
+  // Here we assume that all the models in the chain of models are of type
+  // PhasedArraySpectrumPropagationLossModel that provides the implementation of
+  // this function, i.e. has phased array model of TX and RX as parameters
+  Ptr<SpectrumValue> rxPsd = DoCalcRxPowerSpectralDensity (txPsd, a, b, aPhasedArrayModel, bPhasedArrayModel);
+  if (m_next != 0)
+    {
+      rxPsd = m_next->CalcRxPowerSpectralDensity (rxPsd, a, b, aPhasedArrayModel, bPhasedArrayModel);
+    }
+  return rxPsd;
+}
+
+} // namespace ns3
diff --git a/src/spectrum/model/phased-array-spectrum-propagation-loss-model.h b/src/spectrum/model/phased-array-spectrum-propagation-loss-model.h
new file mode 100644
index 000000000..19c0173b0
--- /dev/null
+++ b/src/spectrum/model/phased-array-spectrum-propagation-loss-model.h
@@ -0,0 +1,120 @@
+/* -*- Mode:C++; c-file-style:"gnu"; indent-tabs-mode:nil; -*- */
+/*
+ * Copyright (c) 2021 CTTC
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License version 2 as
+ * published by the Free Software Foundation;
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ *
+ */
+
+#ifndef PHASED_ARRAY_SPECTRUM_PROPAGATION_LOSS_MODEL_H
+#define PHASED_ARRAY_SPECTRUM_PROPAGATION_LOSS_MODEL_H
+
+
+#include <ns3/object.h>
+#include <ns3/mobility-model.h>
+#include <ns3/spectrum-value.h>
+#include <ns3/phased-array-model.h>
+
+namespace ns3 {
+
+
+/**
+ * \ingroup spectrum
+ *
+ * \brief spectrum-aware propagation loss model that is
+ * compatible with PhasedArrayModel type of ns-3 antenna
+ *
+ * Interface for propagation loss models to be adopted when
+ * transmissions are modeled with a power spectral density by means of
+ * the SpectrumValue class, and when PhasedArrayModel type of atenna
+ * is being used for TX and RX.
+ *
+ */
+class PhasedArraySpectrumPropagationLossModel : public Object
+{
+public:
+  PhasedArraySpectrumPropagationLossModel ();
+  virtual ~PhasedArraySpectrumPropagationLossModel ();
+
+  /**
+   * \brief Get the type ID.
+   * \return the object TypeId
+   */
+  static TypeId GetTypeId ();
+
+  /**
+   * Used to chain various instances of PhasedArraySpectrumPropagationLossModel
+   *
+   * @param next
+   */
+  void SetNext (Ptr<PhasedArraySpectrumPropagationLossModel> next);
+
+  /**
+   * This method is to be called to calculate
+   *
+   * @param txPsd the SpectrumValue representing the power spectral
+   * density of the transmission. Watt units are to be used for radio
+   * communications, and Pascal units for acoustic communications
+   * (e.g., underwater).
+   *
+   * @param a sender mobility
+   * @param b receiver mobility
+   * @param aPhasedArrayModel the instance of the phased antenna array of the sender
+   * @param bPhasedArrayModel the instance of the phased antenna array of the receiver
+   *
+   * @return set of values Vs frequency representing the received
+   * power in the same units used for the txPower parameter.
+   */
+  Ptr<SpectrumValue> CalcRxPowerSpectralDensity (Ptr<const SpectrumValue> txPsd,
+                                                 Ptr<const MobilityModel> a,
+                                                 Ptr<const MobilityModel> b,
+                                                 Ptr<const PhasedArrayModel> aPhasedArrayModel,
+                                                 Ptr<const PhasedArrayModel> bPhasedArrayModel) const;
+
+protected:
+  virtual void DoDispose ();
+
+
+private:
+  /**
+   *
+   * @param txPsd set of values Vs frequency representing the
+   * transmission power. See SpectrumChannel for details.
+   * @param a sender mobility
+   * @param b receiver mobility
+   * @param aPhasedArrayModel the instance of the phased antenna array of the sender
+   * @param bPhasedArrayModel the instance of the phased antenna array of the receiver
+   *
+   * @return set of values Vs frequency representing the received
+   * power in the same units used for the txPower parameter.
+   */
+  virtual Ptr<SpectrumValue> DoCalcRxPowerSpectralDensity (Ptr<const SpectrumValue> txPsd,
+                                                           Ptr<const MobilityModel> a,
+                                                           Ptr<const MobilityModel> b,
+                                                           Ptr<const PhasedArrayModel> aPhasedArrayModel,
+                                                           Ptr<const PhasedArrayModel> bPhasedArrayModel) const = 0;
+
+  Ptr<PhasedArraySpectrumPropagationLossModel> m_next; //!< PhasedArraySpectrumPropagationLossModel chained to this one.
+};
+
+
+
+
+
+
+} // namespace ns3
+
+#endif /* PHASED_ARRAY_SPECTRUM_PROPAGATION_LOSS_MODEL_H */
+
+
diff --git a/src/spectrum/model/single-model-spectrum-channel.cc b/src/spectrum/model/single-model-spectrum-channel.cc
index c572def02..02ed69801 100644
--- a/src/spectrum/model/single-model-spectrum-channel.cc
+++ b/src/spectrum/model/single-model-spectrum-channel.cc
@@ -110,6 +110,19 @@ SingleModelSpectrumChannel::StartTx (Ptr<SpectrumSignalParameters> txParams)
        rxPhyIterator != m_phyList.end ();
        ++rxPhyIterator)
     {
+      Ptr<NetDevice> rxNetDevice = (*rxPhyIterator)->GetDevice ();
+      Ptr<NetDevice> txNetDevice = txParams->txPhy->GetDevice ();
+
+      if (rxNetDevice && txNetDevice)
+        {
+          // we assume that devices are attached to a node
+          if (rxNetDevice->GetNode()->GetId() == txNetDevice->GetNode()->GetId())
+            {
+              NS_LOG_DEBUG ("Skipping the pathloss calculation among different antennas of the same node, not supported yet by any pathloss model in ns-3.");
+              continue;
+            }
+        }
+
       if ((*rxPhyIterator) != txParams->txPhy)
         {
           Time delay  = MicroSeconds (0);
@@ -170,11 +183,10 @@ SingleModelSpectrumChannel::StartTx (Ptr<SpectrumSignalParameters> txParams)
             }
 
 
-          Ptr<NetDevice> netDev = (*rxPhyIterator)->GetDevice ();
-          if (netDev)
+          if (rxNetDevice)
             {
               // the receiver has a NetDevice, so we expect that it is attached to a Node
-              uint32_t dstNode =  netDev->GetNode ()->GetId ();
+              uint32_t dstNode =  rxNetDevice->GetNode ()->GetId ();
               Simulator::ScheduleWithContext (dstNode, delay, &SingleModelSpectrumChannel::StartRx, this, rxParams, *rxPhyIterator);
             }
           else
diff --git a/src/spectrum/model/spectrum-channel.cc b/src/spectrum/model/spectrum-channel.cc
index ea174d309..126b4e200 100644
--- a/src/spectrum/model/spectrum-channel.cc
+++ b/src/spectrum/model/spectrum-channel.cc
@@ -135,6 +135,18 @@ SpectrumChannel::AddSpectrumPropagationLossModel (Ptr<SpectrumPropagationLossMod
   m_spectrumPropagationLoss = loss;
 }
 
+void
+SpectrumChannel::AddPhasedArraySpectrumPropagationLossModel (Ptr<PhasedArraySpectrumPropagationLossModel> loss)
+{
+  NS_LOG_FUNCTION (this << loss);
+  if (m_phasedArraySpectrumPropagationLoss)
+    {
+      loss->SetNext (m_phasedArraySpectrumPropagationLoss);
+    }
+  m_phasedArraySpectrumPropagationLoss = loss;
+}
+
+
 void
 SpectrumChannel::SetPropagationDelayModel (Ptr<PropagationDelayModel> delay)
 {
@@ -149,6 +161,13 @@ SpectrumChannel::GetSpectrumPropagationLossModel (void)
   return m_spectrumPropagationLoss;
 }
 
+Ptr<PhasedArraySpectrumPropagationLossModel>
+SpectrumChannel::GetPhasedArraySpectrumPropagationLossModel (void)
+{
+  NS_LOG_FUNCTION (this);
+  return m_phasedArraySpectrumPropagationLoss;
+}
+
 Ptr<PropagationLossModel>
 SpectrumChannel::GetPropagationLossModel ()
 {
diff --git a/src/spectrum/model/spectrum-channel.h b/src/spectrum/model/spectrum-channel.h
index e314f53e8..c86ea1f8a 100644
--- a/src/spectrum/model/spectrum-channel.h
+++ b/src/spectrum/model/spectrum-channel.h
@@ -27,6 +27,7 @@
 #include <ns3/channel.h>
 #include <ns3/spectrum-signal-parameters.h>
 #include <ns3/spectrum-propagation-loss-model.h>
+#include <ns3/phased-array-spectrum-propagation-loss-model.h>
 #include <ns3/propagation-delay-model.h>
 #include <ns3/propagation-loss-model.h>
 #include <ns3/spectrum-phy.h>
@@ -83,6 +84,13 @@ public:
    */
   void AddSpectrumPropagationLossModel (Ptr<SpectrumPropagationLossModel> loss);
 
+  /**
+   * Add the frequency-dependent propagation loss model
+   * that is compapatible with the phased antenna arrays at the TX and RX
+   * \param loss a pointer to the propagation loss model to be used.
+   */
+  void AddPhasedArraySpectrumPropagationLossModel (Ptr<PhasedArraySpectrumPropagationLossModel> loss);
+
   /**
    * Set the propagation delay model to be used
    * \param delay Ptr to the propagation delay model to be used.
@@ -95,6 +103,13 @@ public:
    */
   Ptr<SpectrumPropagationLossModel> GetSpectrumPropagationLossModel (void);
 
+  /**
+   * Get the frequency-dependent propagation loss model that is
+   * compatible with the phased antenna arrays at TX and RX
+   * \returns a pointer to the propagation loss model.
+   */
+  Ptr<PhasedArraySpectrumPropagationLossModel> GetPhasedArraySpectrumPropagationLossModel (void);
+
   /**
    * Get the propagation loss model.
    * \returns a pointer to the propagation loss model.
@@ -198,6 +213,11 @@ protected:
    */
   Ptr<SpectrumPropagationLossModel> m_spectrumPropagationLoss;
 
+  /**
+   * Frequency-dependent propagation loss model to be used with this channel.
+   */
+  Ptr<PhasedArraySpectrumPropagationLossModel> m_phasedArraySpectrumPropagationLoss;
+
 
 };
 
diff --git a/src/spectrum/model/spectrum-phy.cc b/src/spectrum/model/spectrum-phy.cc
index ae0583d09..6b1aabf0d 100644
--- a/src/spectrum/model/spectrum-phy.cc
+++ b/src/spectrum/model/spectrum-phy.cc
@@ -52,5 +52,24 @@ SpectrumPhy::~SpectrumPhy ()
   NS_LOG_FUNCTION (this);
 }
 
+TypeId
+PhasedArraySpectrumPhy::GetTypeId (void)
+{
+  static TypeId tid = TypeId ("ns3::PhasedArraySpectrumPhy")
+    .SetParent<SpectrumPhy> ()
+    .SetGroupName ("Spectrum")
+  ;
+  return tid;
+}
+
+PhasedArraySpectrumPhy::PhasedArraySpectrumPhy ():SpectrumPhy ()
+{
+  NS_LOG_FUNCTION (this);
+}
+
+PhasedArraySpectrumPhy::~PhasedArraySpectrumPhy ()
+{
+  NS_LOG_FUNCTION (this);
+}
 
 } // namespace
diff --git a/src/spectrum/model/spectrum-phy.h b/src/spectrum/model/spectrum-phy.h
index d02a988ae..4091b2517 100644
--- a/src/spectrum/model/spectrum-phy.h
+++ b/src/spectrum/model/spectrum-phy.h
@@ -31,6 +31,7 @@ class PacketBurst;
 class SpectrumChannel;
 class MobilityModel;
 class AntennaModel;
+class PhasedArrayModel;
 class SpectrumValue;
 class SpectrumModel;
 class NetDevice;
@@ -117,6 +118,49 @@ public:
   virtual void StartRx (Ptr<SpectrumSignalParameters> params) = 0;
 };
 
+/**
+ * \ingroup spectrum
+ *
+ * Abstract base class for Spectrum-aware PHY layers that use
+ * PhasedArrayModel type of antenna for TX/RX
+ *
+ */
+class PhasedArraySpectrumPhy  : public SpectrumPhy
+{
+
+public:
+  PhasedArraySpectrumPhy ();
+  virtual ~PhasedArraySpectrumPhy ();
+
+  /**
+   * \brief Get the type ID.
+   * \return the object TypeId
+   */
+  static TypeId GetTypeId (void);
+
+  /**
+   * Get the PhasedArrayModel used by the NetDevice for TX/RX.
+   *
+   * @return a Ptr to the PhasedArrayModel used by the NetDevice for TX/RX
+   */
+  virtual Ptr<const PhasedArrayModel> GetPhasedArrayModel () const = 0;
+
+private:
+  /**
+   * \brief Copy constructor
+   *
+   * Defined and unimplemented to avoid misuse
+   */
+  PhasedArraySpectrumPhy (PhasedArraySpectrumPhy const &);
+  /**
+   * \brief Copy constructor
+   *
+   * Defined and unimplemented to avoid misuse
+   * \returns
+   */
+  PhasedArraySpectrumPhy& operator= (PhasedArraySpectrumPhy const &);
+};
+
 } // namespace ns3
 
 #endif /* SPECTRUM_PHY_H */
diff --git a/src/spectrum/model/three-gpp-channel-model.cc b/src/spectrum/model/three-gpp-channel-model.cc
index a146dd4c1..b9704963e 100644
--- a/src/spectrum/model/three-gpp-channel-model.cc
+++ b/src/spectrum/model/three-gpp-channel-model.cc
@@ -1064,9 +1064,7 @@ ThreeGppChannelModel::GetChannel (Ptr<const MobilityModel> aMob,
   NS_LOG_FUNCTION (this);
 
   // Compute the channel key. The key is reciprocal, i.e., key (a, b) = key (b, a)
-  uint32_t x1 = std::min (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
-  uint32_t x2 = std::max (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
-  uint32_t channelId = GetKey (x1, x2);
+  uint32_t channelId = GetKey (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
 
   // retrieve the channel condition
   Ptr<const ChannelCondition> condition = m_channelConditionModel->GetChannelCondition (aMob, bMob);
@@ -1647,6 +1645,7 @@ ThreeGppChannelModel::GetNewChannel (Vector locUT, Ptr<const ChannelCondition> c
                       std::tie (rxFieldPatternPhi, rxFieldPatternTheta) = uAntenna->GetElementFieldPattern (Angles (rayAoaRadian[nIndex][mIndex], rayZoaRadian[nIndex][mIndex]));
                       std::tie (txFieldPatternPhi, txFieldPatternTheta) = sAntenna->GetElementFieldPattern (Angles (rayAodRadian[nIndex][mIndex], rayZodRadian[nIndex][mIndex]));
 
+                      NS_ASSERT (4 <= initialPhase.size ());
                       rays += (std::complex<double> (cos (initialPhase[0]), sin (initialPhase[0])) * rxFieldPatternTheta * txFieldPatternTheta +
                                std::complex<double> (cos (initialPhase[1]), sin (initialPhase[1])) * std::sqrt (1 / k) * rxFieldPatternTheta * txFieldPatternPhi +
                                std::complex<double> (cos (initialPhase[2]), sin (initialPhase[2])) * std::sqrt (1 / k) * rxFieldPatternPhi * txFieldPatternTheta +
diff --git a/src/spectrum/model/three-gpp-channel-model.h b/src/spectrum/model/three-gpp-channel-model.h
index c4fb63e9c..462ae755a 100644
--- a/src/spectrum/model/three-gpp-channel-model.h
+++ b/src/spectrum/model/three-gpp-channel-model.h
@@ -256,7 +256,7 @@ private:
    */
   bool ChannelMatrixNeedsUpdate (Ptr<const ThreeGppChannelMatrix> channelMatrix, Ptr<const ChannelCondition> channelCondition) const;
 
-  std::unordered_map<uint32_t, Ptr<ThreeGppChannelMatrix> > m_channelMap; //!< map containing the channel realizations
+  std::unordered_map<uint64_t, Ptr<ThreeGppChannelMatrix> > m_channelMap; //!< map containing the channel realizations
   Time m_updatePeriod; //!< the channel update period
   double m_frequency; //!< the operating frequency
   std::string m_scenario; //!< the 3GPP scenario
diff --git a/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.cc b/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.cc
index 97e7bff45..869289d1c 100644
--- a/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.cc
+++ b/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.cc
@@ -52,7 +52,6 @@ ThreeGppSpectrumPropagationLossModel::~ThreeGppSpectrumPropagationLossModel ()
 void
 ThreeGppSpectrumPropagationLossModel::DoDispose ()
 {
-  m_deviceAntennaMap.clear ();
   m_longTermMap.clear ();
   m_channelModel->Dispose ();
   m_channelModel = nullptr;
@@ -62,7 +61,7 @@ TypeId
 ThreeGppSpectrumPropagationLossModel::GetTypeId (void)
 {
   static TypeId tid = TypeId ("ns3::ThreeGppSpectrumPropagationLossModel")
-    .SetParent<SpectrumPropagationLossModel> ()
+    .SetParent<PhasedArraySpectrumPropagationLossModel> ()
     .SetGroupName ("Spectrum")
     .AddConstructor<ThreeGppSpectrumPropagationLossModel> ()
     .AddAttribute ("ChannelModel",
@@ -94,14 +93,6 @@ ThreeGppSpectrumPropagationLossModel::GetChannelModel () const
   return m_channelModel;
 }
 
-void
-ThreeGppSpectrumPropagationLossModel::AddDevice (Ptr<NetDevice> n, Ptr<const PhasedArrayModel> a)
-{
-  NS_ASSERT_MSG (m_deviceAntennaMap.find (n->GetNode ()->GetId ()) == m_deviceAntennaMap.end (),
-                 "Device is already present in the map");
-  m_deviceAntennaMap.insert (std::make_pair (n->GetNode ()->GetId (), a));
-}
-
 double
 ThreeGppSpectrumPropagationLossModel::GetFrequency () const
 {
@@ -253,9 +244,7 @@ ThreeGppSpectrumPropagationLossModel::GetLongTerm (uint32_t aId, uint32_t bId,
     }
 
   // compute the long term key, the key is unique for each tx-rx pair
-  uint32_t x1 = std::min (aId, bId);
-  uint32_t x2 = std::max (aId, bId);
-  uint32_t longTermId = MatrixBasedChannelModel::GetKey (x1, x2);
+  uint64_t longTermId = MatrixBasedChannelModel::GetKey (aId, bId);
 
   bool update = false; // indicates whether the long term has to be updated
   bool notFound = false; // indicates if the long term has not been computed yet
@@ -302,7 +291,9 @@ ThreeGppSpectrumPropagationLossModel::GetLongTerm (uint32_t aId, uint32_t bId,
 Ptr<SpectrumValue>
 ThreeGppSpectrumPropagationLossModel::DoCalcRxPowerSpectralDensity (Ptr<const SpectrumValue> txPsd,
                                                                     Ptr<const MobilityModel> a,
-                                                                    Ptr<const MobilityModel> b) const
+                                                                    Ptr<const MobilityModel> b,
+                                                                    Ptr<const PhasedArrayModel> aPhasedArrayModel,
+                                                                    Ptr<const PhasedArrayModel> bPhasedArrayModel) const
 {
   NS_LOG_FUNCTION (this);
   uint32_t aId = a->GetObject<Node> ()->GetId (); // id of the node a
@@ -314,20 +305,18 @@ ThreeGppSpectrumPropagationLossModel::DoCalcRxPowerSpectralDensity (Ptr<const Sp
   Ptr<SpectrumValue> rxPsd = Copy<SpectrumValue> (txPsd);
 
   // retrieve the antenna of device a
-  NS_ASSERT_MSG (m_deviceAntennaMap.find (aId) != m_deviceAntennaMap.end (), "Antenna not found for node " << aId);
-  Ptr<const PhasedArrayModel> aAntenna = m_deviceAntennaMap.at (aId);
-  NS_LOG_DEBUG ("a node " << a->GetObject<Node> () << " antenna " << aAntenna);
+  NS_ASSERT_MSG (aPhasedArrayModel, "Antenna not found for node " << aId);
+  NS_LOG_DEBUG ("a node " << a->GetObject<Node> () << " antenna " << aPhasedArrayModel);
 
   // retrieve the antenna of the device b
-  NS_ASSERT_MSG (m_deviceAntennaMap.find (bId) != m_deviceAntennaMap.end (), "Antenna not found for device " << bId);
-  Ptr<const PhasedArrayModel> bAntenna = m_deviceAntennaMap.at (bId);
-  NS_LOG_DEBUG ("b node " << bId << " antenna " << bAntenna);
+  NS_ASSERT_MSG (bPhasedArrayModel, "Antenna not found for device " << bId);
+  NS_LOG_DEBUG ("b node " << bId << " antenna " << bPhasedArrayModel);
 
-  Ptr<const MatrixBasedChannelModel::ChannelMatrix> channelMatrix = m_channelModel->GetChannel (a, b, aAntenna, bAntenna);
+  Ptr<const MatrixBasedChannelModel::ChannelMatrix> channelMatrix = m_channelModel->GetChannel (a, b, aPhasedArrayModel, bPhasedArrayModel);
 
   // get the precoding and combining vectors
-  PhasedArrayModel::ComplexVector aW = aAntenna->GetBeamformingVector ();
-  PhasedArrayModel::ComplexVector bW = bAntenna->GetBeamformingVector ();
+  PhasedArrayModel::ComplexVector aW = aPhasedArrayModel->GetBeamformingVector ();
+  PhasedArrayModel::ComplexVector bW = bPhasedArrayModel->GetBeamformingVector ();
 
   // retrieve the long term component
   PhasedArrayModel::ComplexVector longTerm = GetLongTerm (aId, bId, channelMatrix, aW, bW);
diff --git a/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.h b/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.h
index fb6be13bb..10809a383 100644
--- a/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.h
+++ b/src/spectrum/model/three-gpp-spectrum-propagation-loss-model.h
@@ -22,7 +22,7 @@
 #ifndef THREE_GPP_SPECTRUM_PROPAGATION_LOSS_H
 #define THREE_GPP_SPECTRUM_PROPAGATION_LOSS_H
 
-#include "ns3/spectrum-propagation-loss-model.h"
+#include "ns3/phased-array-spectrum-propagation-loss-model.h"
 #include <complex.h>
 #include <map>
 #include <unordered_map>
@@ -49,7 +49,7 @@ class ChannelCondition;
  * \see PhasedArrayModel
  * \see ChannelCondition
  */
-class ThreeGppSpectrumPropagationLossModel : public SpectrumPropagationLossModel
+class ThreeGppSpectrumPropagationLossModel : public PhasedArraySpectrumPropagationLossModel
 {
 public:
   /**
@@ -82,13 +82,6 @@ public:
    */
   Ptr<MatrixBasedChannelModel> GetChannelModel () const;
 
-  /**
-   * Add a device-antenna pair
-   * \param n a pointer to the NetDevice
-   * \param a a pointer to the associated PhasedArrayModel
-   */
-  void AddDevice (Ptr<NetDevice> n, Ptr<const PhasedArrayModel> a);
-
   /**
    * Sets the value of an attribute belonging to the associated
    * MatrixBasedChannelModel instance
@@ -122,12 +115,15 @@ public:
    * \param txPsd tx PSD
    * \param a first node mobility model
    * \param b second node mobility model
-   *
+   * \param aPhasedArrayModel the antenna array of the first node
+   * \param bPhasedArrayModel the antenna array of the second node
    * \return the received PSD
    */
   Ptr<SpectrumValue> DoCalcRxPowerSpectralDensity (Ptr<const SpectrumValue> txPsd,
                                                    Ptr<const MobilityModel> a,
-                                                   Ptr<const MobilityModel> b) const override;
+                                                   Ptr<const MobilityModel> b,
+                                                   Ptr<const PhasedArrayModel> aPhasedArrayModel,
+                                                   Ptr<const PhasedArrayModel> bPhasedArrayModel) const override;
 
 private:
   /**
@@ -188,8 +184,7 @@ private:
                                           Ptr<const MatrixBasedChannelModel::ChannelMatrix> params,
                                           const Vector &sSpeed, const Vector &uSpeed) const;
 
-  std::unordered_map <uint32_t, Ptr<const PhasedArrayModel> > m_deviceAntennaMap; //!< map containig the <node, antenna> associations
-  mutable std::unordered_map < uint32_t, Ptr<const LongTerm> > m_longTermMap; //!< map containing the long term components
+  mutable std::unordered_map < uint64_t, Ptr<const LongTerm> > m_longTermMap; //!< map containing the long term components
   Ptr<MatrixBasedChannelModel> m_channelModel; //!< the model to generate the channel matrix
 
   // Variable used to compute the additional Doppler contribution for the delayed
diff --git a/src/spectrum/test/three-gpp-channel-test-suite.cc b/src/spectrum/test/three-gpp-channel-test-suite.cc
index 55f05b4e1..39a584184 100644
--- a/src/spectrum/test/three-gpp-channel-test-suite.cc
+++ b/src/spectrum/test/three-gpp-channel-test-suite.cc
@@ -345,6 +345,45 @@ ThreeGppChannelMatrixUpdateTest::DoRun (void)
   Simulator::Destroy ();
 }
 
+/**
+ * \brief A structure that holds the parameters for the function
+ * CheckLongTermUpdate. In this way the problem with the limited
+ * number of parameters of method Schedule is avoided.
+ */
+struct CheckLongTermUpdateParams
+{
+  Ptr<ThreeGppSpectrumPropagationLossModel> lossModel; //!< the ThreeGppSpectrumPropagationLossModel object used to compute the rx PSD
+  Ptr<SpectrumValue> txPsd; //!< the PSD of the tx signal
+  Ptr<MobilityModel> txMob; //!< the mobility model of the tx device
+  Ptr<MobilityModel> rxMob; //!< the mobility model of the rx device
+  Ptr<SpectrumValue> rxPsdOld; //!< the previously received PSD
+  Ptr<PhasedArrayModel> txAntenna; //!< the antenna array of the tx device
+  Ptr<PhasedArrayModel> rxAntenna; //!< the antenna array of the rx device
+
+  /**
+   * \brief Constructor
+   * \param pLossModel the ThreeGppSpectrumPropagationLossModel object used to compute the rx PSD
+   * \param pTxPsd the PSD of the tx signal
+   * \param pTxMob the tx mobility model
+   * \param pRxMob the rx mobility model
+   * \param pRxPsdOld the previously received PSD
+   * \param pTxAntenna the tx antenna array
+   * \param pRxAntenna the rx antenna array
+   */
+  CheckLongTermUpdateParams (Ptr<ThreeGppSpectrumPropagationLossModel> pLossModel, Ptr<SpectrumValue> pTxPsd,
+                             Ptr<MobilityModel> pTxMob, Ptr<MobilityModel> pRxMob, Ptr<SpectrumValue> pRxPsdOld,
+                             Ptr<PhasedArrayModel> pTxAntenna, Ptr<PhasedArrayModel> pRxAntenna)
+  {
+    lossModel = pLossModel;
+    txPsd = pTxPsd;
+    txMob = pTxMob;
+    rxMob = pRxMob;
+    rxPsdOld = pRxPsdOld;
+    txAntenna = pTxAntenna;
+    rxAntenna = pRxAntenna;
+  }
+};
+
 /**
  * Test case for the ThreeGppSpectrumPropagationLossModelTest class.
  * 1) checks if the long term components for the direct and the reverse link
@@ -384,14 +423,9 @@ private:
   /**
    * Test of the long term component is correctly updated when the channel
    * matrix is recomputed
-   * \param lossModel the ThreeGppSpectrumPropagationLossModel object used to
-   *        compute the rx PSD
-   * \param txPsd the PSD of the transmitted signal
-   * \param txMob the mobility model of the tx device
-   * \param rxMob the mobility model of the rx device
-   * \param rxPsdOld the previously received PSD
+   * \param params a structure that contains the set of parameters needed by CheckLongTermUpdate in order to perform calculations
    */
-  void CheckLongTermUpdate (Ptr<ThreeGppSpectrumPropagationLossModel> lossModel, Ptr<SpectrumValue> txPsd, Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, Ptr<SpectrumValue> rxPsdOld);
+  void CheckLongTermUpdate (CheckLongTermUpdateParams &params);
 
   /**
    * Checks if two PSDs are equal
@@ -440,10 +474,10 @@ ThreeGppSpectrumPropagationLossModelTest::ArePsdEqual (Ptr<SpectrumValue> first,
 }
 
 void
-ThreeGppSpectrumPropagationLossModelTest::CheckLongTermUpdate (Ptr<ThreeGppSpectrumPropagationLossModel> lossModel, Ptr<SpectrumValue> txPsd, Ptr<MobilityModel> txMob, Ptr<MobilityModel> rxMob, Ptr<SpectrumValue> rxPsdOld)
+ThreeGppSpectrumPropagationLossModelTest::CheckLongTermUpdate (CheckLongTermUpdateParams &params)
 {
-  Ptr<SpectrumValue> rxPsdNew = lossModel->DoCalcRxPowerSpectralDensity (txPsd, txMob, rxMob);
-  NS_TEST_ASSERT_MSG_EQ (ArePsdEqual (rxPsdOld, rxPsdNew),  false, "The long term is not updated when the channel matrix is recomputed");
+  Ptr<SpectrumValue> rxPsdNew = params.lossModel->DoCalcRxPowerSpectralDensity (params.txPsd, params.txMob, params.rxMob, params.txAntenna, params.rxAntenna);
+  NS_TEST_ASSERT_MSG_EQ (ArePsdEqual (params.rxPsdOld, rxPsdNew),  false, "The long term is not updated when the channel matrix is recomputed");
 }
 
 void
@@ -498,10 +532,6 @@ ThreeGppSpectrumPropagationLossModelTest::DoRun ()
                                                                                     "NumRows", UintegerValue (rxAntennaElements [1]),
                                                                                     "AntennaElement", PointerValue(CreateObject<IsotropicAntennaModel> ()));
   
-  // initialize ThreeGppSpectrumPropagationLossModel
-  lossModel->AddDevice (txDev, txAntenna);
-  lossModel->AddDevice (rxDev, rxAntenna);
-
   // set the beamforming vectors
   DoBeamforming (txDev, txAntenna, rxDev, rxAntenna);
   DoBeamforming (rxDev, rxAntenna, txDev, txAntenna);
@@ -513,10 +543,10 @@ ThreeGppSpectrumPropagationLossModelTest::DoRun ()
   Ptr<SpectrumValue> txPsd =  sf.CreateTxPowerSpectralDensity (txPower, channelNumber);
 
   // compute the rx psd
-  Ptr<SpectrumValue> rxPsdOld = lossModel->DoCalcRxPowerSpectralDensity (txPsd, txMob, rxMob);
+  Ptr<SpectrumValue> rxPsdOld = lossModel->DoCalcRxPowerSpectralDensity (txPsd, txMob, rxMob, txAntenna, rxAntenna);
 
   // 1) check that the rx PSD is equal for both the direct and the reverse channel
-  Ptr<SpectrumValue> rxPsdNew = lossModel->DoCalcRxPowerSpectralDensity (txPsd, rxMob, txMob);
+  Ptr<SpectrumValue> rxPsdNew = lossModel->DoCalcRxPowerSpectralDensity (txPsd, rxMob, txMob, rxAntenna, txAntenna);
   NS_TEST_ASSERT_MSG_EQ (ArePsdEqual (rxPsdOld, rxPsdNew),  true, "The long term for the direct and the reverse channel are different");
 
   // 2) check if the long term is updated when changing the BF vector
@@ -526,7 +556,7 @@ ThreeGppSpectrumPropagationLossModelTest::DoRun ()
   txBfVector [0] = std::complex<double> (0.0, 0.0);
   txAntenna->SetBeamformingVector (txBfVector);
 
-  rxPsdNew = lossModel->DoCalcRxPowerSpectralDensity (txPsd, rxMob, txMob);
+  rxPsdNew = lossModel->DoCalcRxPowerSpectralDensity (txPsd, rxMob, txMob, rxAntenna, txAntenna);
   NS_TEST_ASSERT_MSG_EQ (ArePsdEqual (rxPsdOld, rxPsdNew),  false, "Changing the BF vectors the rx PSD does not change");
 
   // update rxPsdOld
@@ -534,7 +564,7 @@ ThreeGppSpectrumPropagationLossModelTest::DoRun ()
 
   // 3) check if the long term is updated when the channel matrix is recomputed
   Simulator::Schedule (MilliSeconds (101), &ThreeGppSpectrumPropagationLossModelTest::CheckLongTermUpdate,
-                       this, lossModel, txPsd, txMob, rxMob, rxPsdOld);
+                       this, CheckLongTermUpdateParams (lossModel, txPsd, txMob, rxMob, rxPsdOld, txAntenna, rxAntenna));
 
   Simulator::Run ();
   Simulator::Destroy ();
diff --git a/src/spectrum/wscript b/src/spectrum/wscript
index bbf910457..a7346ff24 100644
--- a/src/spectrum/wscript
+++ b/src/spectrum/wscript
@@ -9,6 +9,7 @@ def build(bld):
         'model/spectrum-converter.cc',
         'model/spectrum-signal-parameters.cc',
         'model/spectrum-propagation-loss-model.cc',
+        'model/phased-array-spectrum-propagation-loss-model.cc',
         'model/friis-spectrum-propagation-loss.cc',
         'model/constant-spectrum-propagation-loss.cc',
         'model/spectrum-phy.cc',
@@ -65,6 +66,7 @@ def build(bld):
         'model/spectrum-converter.h',
         'model/spectrum-signal-parameters.h',
         'model/spectrum-propagation-loss-model.h',
+        'model/phased-array-spectrum-propagation-loss-model.h',
         'model/friis-spectrum-propagation-loss.h',
         'model/constant-spectrum-propagation-loss.h',
         'model/spectrum-phy.h',