A Discrete-Event Network Simulator
API
three-gpp-channel-model.cc
Go to the documentation of this file.
1 /* -*- Mode: C++; c-file-style: "gnu"; indent-tabs-mode:nil; -*- */
2 /*
3  * Copyright (c) 2019 SIGNET Lab, Department of Information Engineering,
4  * University of Padova
5  * Copyright (c) 2015, NYU WIRELESS, Tandon School of Engineering,
6  * New York University
7  *
8  * This program is free software; you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License version 2 as
10  * published by the Free Software Foundation;
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software
19  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20  *
21  */
22 
23 #include "three-gpp-channel-model.h"
24 #include "ns3/log.h"
25 #include "ns3/three-gpp-antenna-array-model.h"
26 #include "ns3/node.h"
27 #include "ns3/double.h"
28 #include "ns3/string.h"
29 #include "ns3/integer.h"
30 #include <algorithm>
31 #include <random>
32 #include "ns3/log.h"
33 #include <ns3/simulator.h>
34 #include "ns3/mobility-model.h"
35 #include "ns3/pointer.h"
36 
37 namespace ns3 {
38 
39 NS_LOG_COMPONENT_DEFINE ("ThreeGppChannelModel");
40 
41 NS_OBJECT_ENSURE_REGISTERED (ThreeGppChannelModel);
42 
43 //Table 7.5-3: Ray offset angles within a cluster, given for rms angle spread normalized to 1.
44 static const double offSetAlpha[20] = {
45  0.0447,-0.0447,0.1413,-0.1413,0.2492,-0.2492,0.3715,-0.3715,0.5129,-0.5129,0.6797,-0.6797,0.8844,-0.8844,1.1481,-1.1481,1.5195,-1.5195,2.1551,-2.1551
46 };
47 
48 /*
49  * The cross correlation matrix is constructed according to table 7.5-6.
50  * All the square root matrix is being generated using the Cholesky decomposition
51  * and following the order of [SF,K,DS,ASD,ASA,ZSD,ZSA].
52  * The parameter K is ignored in NLOS.
53  *
54  * The Matlab file to generate the matrices can be found in
55  * https://github.com/nyuwireless-unipd/ns3-mmwave/blob/master/src/mmwave/model/BeamFormingMatrix/SqrtMatrix.m
56  *
57  */
58 static const double sqrtC_RMa_LOS[7][7] = {
59  {1, 0, 0, 0, 0, 0, 0},
60  {0, 1, 0, 0, 0, 0, 0},
61  {-0.5, 0, 0.866025, 0, 0, 0, 0},
62  {0, 0, 0, 1, 0, 0, 0},
63  {0, 0, 0, 0, 1, 0, 0},
64  {0.01, 0, -0.0519615, 0.73, -0.2, 0.651383, 0},
65  {-0.17, -0.02, 0.21362, -0.14, 0.24, 0.142773, 0.909661},
66 };
67 
68 static const double sqrtC_RMa_NLOS[6][6] = {
69  {1, 0, 0, 0, 0, 0},
70  {-0.5, 0.866025, 0, 0, 0, 0},
71  {0.6, -0.11547, 0.791623, 0, 0, 0},
72  {0, 0, 0, 1, 0, 0},
73  {-0.04, -0.138564, 0.540662, -0.18, 0.809003, 0},
74  {-0.25, -0.606218, -0.240013, 0.26, -0.231685, 0.625392},
75 };
76 
77 static const double sqrtC_RMa_O2I[6][6] = {
78  {1, 0, 0, 0, 0, 0},
79  {0, 1, 0, 0, 0, 0},
80  {0, 0, 1, 0, 0, 0},
81  {0, 0, -0.7, 0.714143, 0, 0},
82  {0, 0, 0.66, -0.123225, 0.741091, 0},
83  {0, 0, 0.47, 0.152631, -0.393194, 0.775373},
84 };
85 
86 static const double sqrtC_UMa_LOS[7][7] = {
87  {1, 0, 0, 0, 0, 0, 0},
88  {0, 1, 0, 0, 0, 0, 0},
89  {-0.4, -0.4, 0.824621, 0, 0, 0, 0},
90  {-0.5, 0, 0.242536, 0.83137, 0, 0, 0},
91  {-0.5, -0.2, 0.630593, -0.484671, 0.278293, 0, 0},
92  {0, 0, -0.242536, 0.672172, 0.642214, 0.27735, 0},
93  {-0.8, 0, -0.388057, -0.367926, 0.238537, -3.58949e-15, 0.130931},
94 };
95 
96 
97 static const double sqrtC_UMa_NLOS[6][6] = {
98  {1, 0, 0, 0, 0, 0},
99  {-0.4, 0.916515, 0, 0, 0, 0},
100  {-0.6, 0.174574, 0.78072, 0, 0, 0},
101  {0, 0.654654, 0.365963, 0.661438, 0, 0},
102  {0, -0.545545, 0.762422, 0.118114, 0.327327, 0},
103  {-0.4, -0.174574, -0.396459, 0.392138, 0.49099, 0.507445},
104 };
105 
106 static const double sqrtC_UMa_O2I[6][6] = {
107  {1, 0, 0, 0, 0, 0},
108  {-0.5, 0.866025, 0, 0, 0, 0},
109  {0.2, 0.57735, 0.791623, 0, 0, 0},
110  {0, 0.46188, -0.336861, 0.820482, 0, 0},
111  {0, -0.69282, 0.252646, 0.493742, 0.460857, 0},
112  {0, -0.23094, 0.16843, 0.808554, -0.220827, 0.464515},
113 
114 };
115 
116 static const double sqrtC_UMi_LOS[7][7] = {
117  {1, 0, 0, 0, 0, 0, 0},
118  {0.5, 0.866025, 0, 0, 0, 0, 0},
119  {-0.4, -0.57735, 0.711805, 0, 0, 0, 0},
120  {-0.5, 0.057735, 0.468293, 0.726201, 0, 0, 0},
121  {-0.4, -0.11547, 0.805464, -0.23482, 0.350363, 0, 0},
122  {0, 0, 0, 0.688514, 0.461454, 0.559471, 0},
123  {0, 0, 0.280976, 0.231921, -0.490509, 0.11916, 0.782603},
124 };
125 
126 static const double sqrtC_UMi_NLOS[6][6] = {
127  {1, 0, 0, 0, 0, 0},
128  {-0.7, 0.714143, 0, 0, 0, 0},
129  {0, 0, 1, 0, 0, 0},
130  {-0.4, 0.168034, 0, 0.90098, 0, 0},
131  {0, -0.70014, 0.5, 0.130577, 0.4927, 0},
132  {0, 0, 0.5, 0.221981, -0.566238, 0.616522},
133 };
134 
135 static const double sqrtC_UMi_O2I[6][6] = {
136  {1, 0, 0, 0, 0, 0},
137  {-0.5, 0.866025, 0, 0, 0, 0},
138  {0.2, 0.57735, 0.791623, 0, 0, 0},
139  {0, 0.46188, -0.336861, 0.820482, 0, 0},
140  {0, -0.69282, 0.252646, 0.493742, 0.460857, 0},
141  {0, -0.23094, 0.16843, 0.808554, -0.220827, 0.464515},
142 };
143 
144 static const double sqrtC_office_LOS[7][7] = {
145  {1, 0, 0, 0, 0, 0, 0},
146  {0.5, 0.866025, 0, 0, 0, 0, 0},
147  {-0.8, -0.11547, 0.588784, 0, 0, 0, 0},
148  {-0.4, 0.23094, 0.520847, 0.717903, 0, 0, 0},
149  {-0.5, 0.288675, 0.73598, -0.348236, 0.0610847, 0, 0},
150  {0.2, -0.11547, 0.418943, 0.541106, 0.219905, 0.655744, 0},
151  {0.3, -0.057735, 0.73598, -0.348236, 0.0610847, -0.304997, 0.383375},
152 };
153 
154 static const double sqrtC_office_NLOS[6][6] = {
155  {1, 0, 0, 0, 0, 0},
156  {-0.5, 0.866025, 0, 0, 0, 0},
157  {0, 0.46188, 0.886942, 0, 0, 0},
158  {-0.4, -0.23094, 0.120263, 0.878751, 0, 0},
159  {0, -0.311769, 0.55697, -0.249198, 0.728344, 0},
160  {0, -0.069282, 0.295397, 0.430696, 0.468462, 0.709214},
161 };
162 
163 ThreeGppChannelModel::ThreeGppChannelModel ()
164 {
165  NS_LOG_FUNCTION (this);
166  m_uniformRv = CreateObject<UniformRandomVariable> ();
167  m_uniformRvShuffle = CreateObject<UniformRandomVariable> ();
168 
169  m_normalRv = CreateObject<NormalRandomVariable> ();
170  m_normalRv->SetAttribute ("Mean", DoubleValue (0.0));
171  m_normalRv->SetAttribute ("Variance", DoubleValue (1.0));
172 }
173 
174 ThreeGppChannelModel::~ThreeGppChannelModel ()
175 {
176  NS_LOG_FUNCTION (this);
177 }
178 
179 void
180 ThreeGppChannelModel::DoDispose ()
181 {
182  m_channelMap.clear ();
183  m_channelConditionModel->Dispose ();
184  m_channelConditionModel = nullptr;
185 }
186 
187 TypeId
188 ThreeGppChannelModel::GetTypeId (void)
189 {
190  static TypeId tid = TypeId ("ns3::ThreeGppChannelModel")
191  .SetParent<Object> ()
192  .SetGroupName ("Spectrum")
193  .SetParent<MatrixBasedChannelModel> ()
194  .AddConstructor<ThreeGppChannelModel> ()
195  .AddAttribute ("Frequency",
196  "The operating Frequency in Hz",
197  DoubleValue (500.0e6),
198  MakeDoubleAccessor (&ThreeGppChannelModel::SetFrequency,
199  &ThreeGppChannelModel::GetFrequency),
200  MakeDoubleChecker<double> ())
201  .AddAttribute ("Scenario",
202  "The 3GPP scenario (RMa, UMa, UMi-StreetCanyon, InH-OfficeOpen, InH-OfficeMixed)",
203  StringValue ("UMa"),
204  MakeStringAccessor (&ThreeGppChannelModel::SetScenario,
205  &ThreeGppChannelModel::GetScenario),
206  MakeStringChecker ())
207  .AddAttribute ("ChannelConditionModel",
208  "Pointer to the channel condition model",
209  PointerValue (),
210  MakePointerAccessor (&ThreeGppChannelModel::SetChannelConditionModel,
211  &ThreeGppChannelModel::GetChannelConditionModel),
212  MakePointerChecker<ChannelConditionModel> ())
213  .AddAttribute ("UpdatePeriod",
214  "Specify the channel coherence time",
215  TimeValue (MilliSeconds (0)),
216  MakeTimeAccessor (&ThreeGppChannelModel::m_updatePeriod),
217  MakeTimeChecker ())
218  // attributes for the blockage model
219  .AddAttribute ("Blockage",
220  "Enable blockage model A (sec 7.6.4.1)",
221  BooleanValue (false),
222  MakeBooleanAccessor (&ThreeGppChannelModel::m_blockage),
223  MakeBooleanChecker ())
224  .AddAttribute ("NumNonselfBlocking",
225  "number of non-self-blocking regions",
226  IntegerValue (4),
227  MakeIntegerAccessor (&ThreeGppChannelModel::m_numNonSelfBlocking),
228  MakeIntegerChecker<uint16_t> ())
229  .AddAttribute ("PortraitMode",
230  "true for portrait mode, false for landscape mode",
231  BooleanValue (true),
232  MakeBooleanAccessor (&ThreeGppChannelModel::m_portraitMode),
233  MakeBooleanChecker ())
234  .AddAttribute ("BlockerSpeed",
235  "The speed of moving blockers, the unit is m/s",
236  DoubleValue (1),
237  MakeDoubleAccessor (&ThreeGppChannelModel::m_blockerSpeed),
238  MakeDoubleChecker<double> ())
239  ;
240  return tid;
241 }
242 
243 void
244 ThreeGppChannelModel::SetChannelConditionModel (Ptr<ChannelConditionModel> model)
245 {
246  NS_LOG_FUNCTION (this);
247  m_channelConditionModel = model;
248 }
249 
250 Ptr<ChannelConditionModel>
251 ThreeGppChannelModel::GetChannelConditionModel () const
252 {
253  NS_LOG_FUNCTION (this);
254  return m_channelConditionModel;
255 }
256 
257 void
258 ThreeGppChannelModel::SetFrequency (double f)
259 {
260  NS_LOG_FUNCTION (this);
261  NS_ASSERT_MSG (f >= 500.0e6 && f <= 100.0e9, "Frequency should be between 0.5 and 100 GHz but is " << f);
262  m_frequency = f;
263 }
264 
265 double
266 ThreeGppChannelModel::GetFrequency () const
267 {
268  NS_LOG_FUNCTION (this);
269  return m_frequency;
270 }
271 
272 void
273 ThreeGppChannelModel::SetScenario (const std::string &scenario)
274 {
275  NS_LOG_FUNCTION (this);
276  NS_ASSERT_MSG (scenario == "RMa" || scenario == "UMa" || scenario == "UMi-StreetCanyon"
277  || scenario == "InH-OfficeOpen" || scenario == "InH-OfficeMixed"
278  || scenario == "V2V-Urban" || scenario == "V2V-Highway",
279  "Unknown scenario, choose between RMa, UMa, UMi-StreetCanyon,"
280  "InH-OfficeOpen, InH-OfficeMixed, V2V-Urban or V2V-Highway");
281  m_scenario = scenario;
282 }
283 
284 std::string
285 ThreeGppChannelModel::GetScenario () const
286 {
287  NS_LOG_FUNCTION (this);
288  return m_scenario;
289 }
290 
291 Ptr<const ThreeGppChannelModel::ParamsTable>
292 ThreeGppChannelModel::GetThreeGppTable (Ptr<const ChannelCondition> channelCondition, double hBS, double hUT, double distance2D) const
293 {
294  NS_LOG_FUNCTION (this);
295 
296  double fcGHz = m_frequency / 1e9;
297  Ptr<ParamsTable> table3gpp = Create<ParamsTable> ();
298  // table3gpp includes the following parameters:
299  // numOfCluster, raysPerCluster, uLgDS, sigLgDS, uLgASD, sigLgASD,
300  // uLgASA, sigLgASA, uLgZSA, sigLgZSA, uLgZSD, sigLgZSD, offsetZOD,
301  // cDS, cASD, cASA, cZSA, uK, sigK, rTau, uXpr, sigXpr, shadowingStd
302 
303  bool los = channelCondition->IsLos ();
304  bool o2i = channelCondition->IsO2i ();
305 
306  // In NLOS case, parameter uK and sigK are not used and they are set to 0
307  if (m_scenario == "RMa")
308  {
309  if (los && !o2i)
310  {
311  // 3GPP mentioned that 3.91 ns should be used when the Cluster DS (cDS)
312  // entry is N/A.
313  table3gpp->m_numOfCluster = 11;
314  table3gpp->m_raysPerCluster = 20;
315  table3gpp->m_uLgDS = -7.49;
316  table3gpp->m_sigLgDS = 0.55;
317  table3gpp->m_uLgASD = 0.90;
318  table3gpp->m_sigLgASD = 0.38;
319  table3gpp->m_uLgASA = 1.52;
320  table3gpp->m_sigLgASA = 0.24;
321  table3gpp->m_uLgZSA = 0.47;
322  table3gpp->m_sigLgZSA = 0.40;
323  table3gpp->m_uLgZSD = 0.34;
324  table3gpp->m_sigLgZSD = std::max (-1.0, -0.17 * (distance2D / 1000) - 0.01 * (hUT - 1.5) + 0.22);
325  table3gpp->m_offsetZOD = 0;
326  table3gpp->m_cDS = 3.91e-9;
327  table3gpp->m_cASD = 2;
328  table3gpp->m_cASA = 3;
329  table3gpp->m_cZSA = 3;
330  table3gpp->m_uK = 7;
331  table3gpp->m_sigK = 4;
332  table3gpp->m_rTau = 3.8;
333  table3gpp->m_uXpr = 12;
334  table3gpp->m_sigXpr = 4;
335  table3gpp->m_perClusterShadowingStd = 3;
336 
337  for (uint8_t row = 0; row < 7; row++)
338  {
339  for (uint8_t column = 0; column < 7; column++)
340  {
341  table3gpp->m_sqrtC[row][column] = sqrtC_RMa_LOS[row][column];
342  }
343  }
344  }
345  else if (!los && !o2i)
346  {
347  table3gpp->m_numOfCluster = 10;
348  table3gpp->m_raysPerCluster = 20;
349  table3gpp->m_uLgDS = -7.43;
350  table3gpp->m_sigLgDS = 0.48;
351  table3gpp->m_uLgASD = 0.95;
352  table3gpp->m_sigLgASD = 0.45;
353  table3gpp->m_uLgASA = 1.52;
354  table3gpp->m_sigLgASA = 0.13;
355  table3gpp->m_uLgZSA = 0.58,
356  table3gpp->m_sigLgZSA = 0.37;
357  table3gpp->m_uLgZSD = std::max (-1.0, -0.19 * (distance2D / 1000) - 0.01 * (hUT - 1.5) + 0.28);
358  table3gpp->m_sigLgZSD = 0.30;
359  table3gpp->m_offsetZOD = atan ((35 - 3.5) / distance2D) - atan ((35 - 1.5) / distance2D);
360  table3gpp->m_cDS = 3.91e-9;
361  table3gpp->m_cASD = 2;
362  table3gpp->m_cASA = 3;
363  table3gpp->m_cZSA = 3;
364  table3gpp->m_uK = 0;
365  table3gpp->m_sigK = 0;
366  table3gpp->m_rTau = 1.7;
367  table3gpp->m_uXpr = 7;
368  table3gpp->m_sigXpr = 3;
369  table3gpp->m_perClusterShadowingStd = 3;
370 
371  for (uint8_t row = 0; row < 6; row++)
372  {
373  for (uint8_t column = 0; column < 6; column++)
374  {
375  table3gpp->m_sqrtC[row][column] = sqrtC_RMa_NLOS[row][column];
376  }
377  }
378  }
379  else // o2i
380  {
381  table3gpp->m_numOfCluster = 10;
382  table3gpp->m_raysPerCluster = 20;
383  table3gpp->m_uLgDS = -7.47;
384  table3gpp->m_sigLgDS = 0.24;
385  table3gpp->m_uLgASD = 0.67;
386  table3gpp->m_sigLgASD = 0.18;
387  table3gpp->m_uLgASA = 1.66;
388  table3gpp->m_sigLgASA = 0.21;
389  table3gpp->m_uLgZSA = 0.93,
390  table3gpp->m_sigLgZSA = 0.22;
391  table3gpp->m_uLgZSD = std::max (-1.0, -0.19 * (distance2D / 1000) - 0.01 * (hUT - 1.5) + 0.28);
392  table3gpp->m_sigLgZSD = 0.30;
393  table3gpp->m_offsetZOD = atan ((35 - 3.5) / distance2D) - atan ((35 - 1.5) / distance2D);
394  table3gpp->m_cDS = 3.91e-9;
395  table3gpp->m_cASD = 2;
396  table3gpp->m_cASA = 3;
397  table3gpp->m_cZSA = 3;
398  table3gpp->m_uK = 0;
399  table3gpp->m_sigK = 0;
400  table3gpp->m_rTau = 1.7;
401  table3gpp->m_uXpr = 7;
402  table3gpp->m_sigXpr = 3;
403  table3gpp->m_perClusterShadowingStd = 3;
404 
405  for (uint8_t row = 0; row < 6; row++)
406  {
407  for (uint8_t column = 0; column < 6; column++)
408  {
409  table3gpp->m_sqrtC[row][column] = sqrtC_RMa_O2I[row][column];
410  }
411  }
412  }
413  }
414  else if (m_scenario == "UMa")
415  {
416  if (los && !o2i)
417  {
418  table3gpp->m_numOfCluster = 12;
419  table3gpp->m_raysPerCluster = 20;
420  table3gpp->m_uLgDS = -6.955 - 0.0963 * log10 (fcGHz);
421  table3gpp->m_sigLgDS = 0.66;
422  table3gpp->m_uLgASD = 1.06 + 0.1114 * log10 (fcGHz);
423  table3gpp->m_sigLgASD = 0.28;
424  table3gpp->m_uLgASA = 1.81;
425  table3gpp->m_sigLgASA = 0.20;
426  table3gpp->m_uLgZSA = 0.95;
427  table3gpp->m_sigLgZSA = 0.16;
428  table3gpp->m_uLgZSD = std::max (-0.5, -2.1 * distance2D / 1000 - 0.01 * (hUT - 1.5) + 0.75);
429  table3gpp->m_sigLgZSD = 0.40;
430  table3gpp->m_offsetZOD = 0;
431  table3gpp->m_cDS = std::max (0.25, -3.4084 * log10 (fcGHz) + 6.5622) * 1e-9;
432  table3gpp->m_cASD = 5;
433  table3gpp->m_cASA = 11;
434  table3gpp->m_cZSA = 7;
435  table3gpp->m_uK = 9;
436  table3gpp->m_sigK = 3.5;
437  table3gpp->m_rTau = 2.5;
438  table3gpp->m_uXpr = 8;
439  table3gpp->m_sigXpr = 4;
440  table3gpp->m_perClusterShadowingStd = 3;
441 
442  for (uint8_t row = 0; row < 7; row++)
443  {
444  for (uint8_t column = 0; column < 7; column++)
445  {
446  table3gpp->m_sqrtC[row][column] = sqrtC_UMa_LOS[row][column];
447  }
448  }
449  }
450  else
451  {
452  double uLgZSD = std::max (-0.5, -2.1 * distance2D / 1000 - 0.01 * (hUT - 1.5) + 0.9);
453 
454  double afc = 0.208 * log10 (fcGHz) - 0.782;
455  double bfc = 25;
456  double cfc = -0.13 * log10 (fcGHz) + 2.03;
457  double efc = 7.66 * log10 (fcGHz) - 5.96;
458 
459  double offsetZOD = efc - std::pow (10, afc * log10 (std::max (bfc,distance2D)) + cfc);
460 
461  if (!los && !o2i)
462  {
463  table3gpp->m_numOfCluster = 20;
464  table3gpp->m_raysPerCluster = 20;
465  table3gpp->m_uLgDS = -6.28 - 0.204 * log10 (fcGHz);
466  table3gpp->m_sigLgDS = 0.39;
467  table3gpp->m_uLgASD = 1.5 - 0.1144 * log10 (fcGHz);
468  table3gpp->m_sigLgASD = 0.28;
469  table3gpp->m_uLgASA = 2.08 - 0.27 * log10 (fcGHz);
470  table3gpp->m_sigLgASA = 0.11;
471  table3gpp->m_uLgZSA = -0.3236 * log10 (fcGHz) + 1.512;
472  table3gpp->m_sigLgZSA = 0.16;
473  table3gpp->m_uLgZSD = uLgZSD;
474  table3gpp->m_sigLgZSD = 0.49;
475  table3gpp->m_offsetZOD = offsetZOD;
476  table3gpp->m_cDS = std::max (0.25, -3.4084 * log10 (fcGHz) + 6.5622) * 1e-9;
477  table3gpp->m_cASD = 2;
478  table3gpp->m_cASA = 15;
479  table3gpp->m_cZSA = 7;
480  table3gpp->m_uK = 0;
481  table3gpp->m_sigK = 0;
482  table3gpp->m_rTau = 2.3;
483  table3gpp->m_uXpr = 7;
484  table3gpp->m_sigXpr = 3;
485  table3gpp->m_perClusterShadowingStd = 3;
486 
487  for (uint8_t row = 0; row < 6; row++)
488  {
489  for (uint8_t column = 0; column < 6; column++)
490  {
491  table3gpp->m_sqrtC[row][column] = sqrtC_UMa_NLOS[row][column];
492  }
493  }
494  }
495  else //(o2i)
496  {
497  table3gpp->m_numOfCluster = 12;
498  table3gpp->m_raysPerCluster = 20;
499  table3gpp->m_uLgDS = -6.62;
500  table3gpp->m_sigLgDS = 0.32;
501  table3gpp->m_uLgASD = 1.25;
502  table3gpp->m_sigLgASD = 0.42;
503  table3gpp->m_uLgASA = 1.76;
504  table3gpp->m_sigLgASA = 0.16;
505  table3gpp->m_uLgZSA = 1.01;
506  table3gpp->m_sigLgZSA = 0.43;
507  table3gpp->m_uLgZSD = uLgZSD;
508  table3gpp->m_sigLgZSD = 0.49;
509  table3gpp->m_offsetZOD = offsetZOD;
510  table3gpp->m_cDS = 11e-9;
511  table3gpp->m_cASD = 5;
512  table3gpp->m_cASA = 8;
513  table3gpp->m_cZSA = 3;
514  table3gpp->m_uK = 0;
515  table3gpp->m_sigK = 0;
516  table3gpp->m_rTau = 2.2;
517  table3gpp->m_uXpr = 9;
518  table3gpp->m_sigXpr = 5;
519  table3gpp->m_perClusterShadowingStd = 4;
520 
521  for (uint8_t row = 0; row < 6; row++)
522  {
523  for (uint8_t column = 0; column < 6; column++)
524  {
525  table3gpp->m_sqrtC[row][column] = sqrtC_UMa_O2I[row][column];
526  }
527  }
528 
529  }
530 
531  }
532 
533  }
534  else if (m_scenario == "UMi-StreetCanyon")
535  {
536  if (los && !o2i)
537  {
538  table3gpp->m_numOfCluster = 12;
539  table3gpp->m_raysPerCluster = 20;
540  table3gpp->m_uLgDS = -0.24 * log10 (1 + fcGHz) - 7.14;
541  table3gpp->m_sigLgDS = 0.38;
542  table3gpp->m_uLgASD = -0.05 * log10 (1 + fcGHz) + 1.21;
543  table3gpp->m_sigLgASD = 0.41;
544  table3gpp->m_uLgASA = -0.08 * log10 (1 + fcGHz) + 1.73;
545  table3gpp->m_sigLgASA = 0.014 * log10 (1 + fcGHz) + 0.28;
546  table3gpp->m_uLgZSA = -0.1 * log10 (1 + fcGHz) + 0.73;
547  table3gpp->m_sigLgZSA = -0.04 * log10 (1 + fcGHz) + 0.34;
548  table3gpp->m_uLgZSD = std::max (-0.21, -14.8 * distance2D / 1000 + 0.01 * std::abs (hUT - hBS) + 0.83);
549  table3gpp->m_sigLgZSD = 0.35;
550  table3gpp->m_offsetZOD = 0;
551  table3gpp->m_cDS = 5e-9;
552  table3gpp->m_cASD = 3;
553  table3gpp->m_cASA = 17;
554  table3gpp->m_cZSA = 7;
555  table3gpp->m_uK = 9;
556  table3gpp->m_sigK = 5;
557  table3gpp->m_rTau = 3;
558  table3gpp->m_uXpr = 9;
559  table3gpp->m_sigXpr = 3;
560  table3gpp->m_perClusterShadowingStd = 3;
561 
562  for (uint8_t row = 0; row < 7; row++)
563  {
564  for (uint8_t column = 0; column < 7; column++)
565  {
566  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
567  }
568  }
569  }
570  else
571  {
572  double uLgZSD = std::max (-0.5, -3.1 * distance2D / 1000 + 0.01 * std::max (hUT - hBS,0.0) + 0.2);
573  double offsetZOD = -1 * std::pow (10, -1.5 * log10 (std::max (10.0, distance2D)) + 3.3);
574  if (!los && !o2i)
575  {
576  table3gpp->m_numOfCluster = 19;
577  table3gpp->m_raysPerCluster = 20;
578  table3gpp->m_uLgDS = -0.24 * log10 (1 + fcGHz) - 6.83;
579  table3gpp->m_sigLgDS = 0.16 * log10 (1 + fcGHz) + 0.28;
580  table3gpp->m_uLgASD = -0.23 * log10 (1 + fcGHz) + 1.53;
581  table3gpp->m_sigLgASD = 0.11 * log10 (1 + fcGHz) + 0.33;
582  table3gpp->m_uLgASA = -0.08 * log10 (1 + fcGHz) + 1.81;
583  table3gpp->m_sigLgASA = 0.05 * log10 (1 + fcGHz) + 0.3;
584  table3gpp->m_uLgZSA = -0.04 * log10 (1 + fcGHz) + 0.92;
585  table3gpp->m_sigLgZSA = -0.07 * log10 (1 + fcGHz) + 0.41;
586  table3gpp->m_uLgZSD = uLgZSD;
587  table3gpp->m_sigLgZSD = 0.35;
588  table3gpp->m_offsetZOD = offsetZOD;
589  table3gpp->m_cDS = 11e-9;
590  table3gpp->m_cASD = 10;
591  table3gpp->m_cASA = 22;
592  table3gpp->m_cZSA = 7;
593  table3gpp->m_uK = 0;
594  table3gpp->m_sigK = 0;
595  table3gpp->m_rTau = 2.1;
596  table3gpp->m_uXpr = 8;
597  table3gpp->m_sigXpr = 3;
598  table3gpp->m_perClusterShadowingStd = 3;
599 
600  for (uint8_t row = 0; row < 6; row++)
601  {
602  for (uint8_t column = 0; column < 6; column++)
603  {
604  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_NLOS[row][column];
605  }
606  }
607  }
608  else //(o2i)
609  {
610  table3gpp->m_numOfCluster = 12;
611  table3gpp->m_raysPerCluster = 20;
612  table3gpp->m_uLgDS = -6.62;
613  table3gpp->m_sigLgDS = 0.32;
614  table3gpp->m_uLgASD = 1.25;
615  table3gpp->m_sigLgASD = 0.42;
616  table3gpp->m_uLgASA = 1.76;
617  table3gpp->m_sigLgASA = 0.16;
618  table3gpp->m_uLgZSA = 1.01;
619  table3gpp->m_sigLgZSA = 0.43;
620  table3gpp->m_uLgZSD = uLgZSD;
621  table3gpp->m_sigLgZSD = 0.35;
622  table3gpp->m_offsetZOD = offsetZOD;
623  table3gpp->m_cDS = 11e-9;
624  table3gpp->m_cASD = 5;
625  table3gpp->m_cASA = 8;
626  table3gpp->m_cZSA = 3;
627  table3gpp->m_uK = 0;
628  table3gpp->m_sigK = 0;
629  table3gpp->m_rTau = 2.2;
630  table3gpp->m_uXpr = 9;
631  table3gpp->m_sigXpr = 5;
632  table3gpp->m_perClusterShadowingStd = 4;
633 
634  for (uint8_t row = 0; row < 6; row++)
635  {
636  for (uint8_t column = 0; column < 6; column++)
637  {
638  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_O2I[row][column];
639  }
640  }
641  }
642  }
643  }
644  else if (m_scenario == "InH-OfficeMixed"||m_scenario == "InH-OfficeOpen")
645  {
646  NS_ASSERT_MSG (!o2i, "The indoor scenario does out support outdoor to indoor");
647  if (los)
648  {
649  table3gpp->m_numOfCluster = 15;
650  table3gpp->m_raysPerCluster = 20;
651  table3gpp->m_uLgDS = -0.01 * log10 (1 + fcGHz) - 7.692;
652  table3gpp->m_sigLgDS = 0.18;
653  table3gpp->m_uLgASD = 1.60;
654  table3gpp->m_sigLgASD = 0.18;
655  table3gpp->m_uLgASA = -0.19 * log10 (1 + fcGHz) + 1.781;
656  table3gpp->m_sigLgASA = 0.12 * log10 (1 + fcGHz) + 0.119;
657  table3gpp->m_uLgZSA = -0.26 * log10 (1 + fcGHz) + 1.44;
658  table3gpp->m_sigLgZSA = -0.04 * log10 (1 + fcGHz) + 0.264;
659  table3gpp->m_uLgZSD = -1.43 * log10 (1 + fcGHz) + 2.228;
660  table3gpp->m_sigLgZSD = 0.13 * log10 (1 + fcGHz) + 0.30;
661  table3gpp->m_offsetZOD = 0;
662  table3gpp->m_cDS = 3.91e-9;
663  table3gpp->m_cASD = 5;
664  table3gpp->m_cASA = 8;
665  table3gpp->m_cZSA = 9;
666  table3gpp->m_uK = 7;
667  table3gpp->m_sigK = 4;
668  table3gpp->m_rTau = 3.6;
669  table3gpp->m_uXpr = 11;
670  table3gpp->m_sigXpr = 4;
671  table3gpp->m_perClusterShadowingStd = 6;
672 
673  for (uint8_t row = 0; row < 7; row++)
674  {
675  for (uint8_t column = 0; column < 7; column++)
676  {
677  table3gpp->m_sqrtC[row][column] = sqrtC_office_LOS[row][column];
678  }
679  }
680  }
681  else
682  {
683  table3gpp->m_numOfCluster = 19;
684  table3gpp->m_raysPerCluster = 20;
685  table3gpp->m_uLgDS = -0.28 * log10 (1 + fcGHz) - 7.173;
686  table3gpp->m_sigLgDS = 0.1 * log10 (1 + fcGHz) + 0.055;
687  table3gpp->m_uLgASD = 1.62;
688  table3gpp->m_sigLgASD = 0.25;
689  table3gpp->m_uLgASA = -0.11 * log10 (1 + fcGHz) + 1.863;
690  table3gpp->m_sigLgASA = 0.12 * log10 (1 + fcGHz) + 0.059;
691  table3gpp->m_uLgZSA = -0.15 * log10 (1 + fcGHz) + 1.387;
692  table3gpp->m_sigLgZSA = -0.09 * log10 (1 + fcGHz) + 0.746;
693  table3gpp->m_uLgZSD = 1.08;
694  table3gpp->m_sigLgZSD = 0.36;
695  table3gpp->m_offsetZOD = 0;
696  table3gpp->m_cDS = 3.91e-9;
697  table3gpp->m_cASD = 5;
698  table3gpp->m_cASA = 11;
699  table3gpp->m_cZSA = 9;
700  table3gpp->m_uK = 0;
701  table3gpp->m_sigK = 0;
702  table3gpp->m_rTau = 3;
703  table3gpp->m_uXpr = 10;
704  table3gpp->m_sigXpr = 4;
705  table3gpp->m_perClusterShadowingStd = 3;
706 
707  for (uint8_t row = 0; row < 6; row++)
708  {
709  for (uint8_t column = 0; column < 6; column++)
710  {
711  table3gpp->m_sqrtC[row][column] = sqrtC_office_NLOS[row][column];
712  }
713  }
714  }
715  }
716  else if (m_scenario == "V2V-Urban")
717  {
718  if (channelCondition->IsLos ())
719  {
720  // 3GPP mentioned that 3.91 ns should be used when the Cluster DS (cDS)
721  // entry is N/A.
722  table3gpp->m_numOfCluster = 12;
723  table3gpp->m_raysPerCluster = 20;
724  table3gpp->m_uLgDS = -0.2 * log10 (1 + fcGHz) - 7.5;
725  table3gpp->m_sigLgDS = 0.1;
726  table3gpp->m_uLgASD = -0.1 * log10 (1 + fcGHz) + 1.6;
727  table3gpp->m_sigLgASD = 0.1;
728  table3gpp->m_uLgASA = -0.1 * log10 (1 + fcGHz) + 1.6;
729  table3gpp->m_sigLgASA = 0.1;
730  table3gpp->m_uLgZSA = -0.1 * log10 (1 + fcGHz) + 0.73,
731  table3gpp->m_sigLgZSA = -0.04 * log10 (1 + fcGHz) + 0.34,
732  table3gpp->m_uLgZSD = -0.1 * log10 (1 + fcGHz) + 0.73;
733  table3gpp->m_sigLgZSD = -0.04 * log10 (1 + fcGHz) + 0.34;
734  table3gpp->m_offsetZOD = 0;
735  table3gpp->m_cDS = 5;
736  table3gpp->m_cASD = 17;
737  table3gpp->m_cASA = 17;
738  table3gpp->m_cZSA = 7;
739  table3gpp->m_uK = 3.48;
740  table3gpp->m_sigK = 2;
741  table3gpp->m_rTau = 3;
742  table3gpp->m_uXpr = 9;
743  table3gpp->m_sigXpr = 3;
744  table3gpp->m_perClusterShadowingStd = 4;
745 
746  for (uint8_t row = 0; row < 7; row++)
747  {
748  for (uint8_t column = 0; column < 7; column++)
749  {
750  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
751  }
752  }
753  }
754  else if (channelCondition->IsNlos ())
755  {
756  table3gpp->m_numOfCluster = 19;
757  table3gpp->m_raysPerCluster = 20;
758  table3gpp->m_uLgDS = -0.3 * log10 (1 + fcGHz) - 7;
759  table3gpp->m_sigLgDS = 0.28;
760  table3gpp->m_uLgASD = -0.08 * log10 (1 + fcGHz) + 1.81;
761  table3gpp->m_sigLgASD = 0.05 * log10 (1 + fcGHz) + 0.3;
762  table3gpp->m_uLgASA = -0.08 * log10 (1 + fcGHz) + 1.81;
763  table3gpp->m_sigLgASA = 0.05 * log10 (1 + fcGHz) + 0.3;
764  table3gpp->m_uLgZSA = -0.04 * log10 (1 + fcGHz) + 0.92,
765  table3gpp->m_sigLgZSA = -0.07 * log10 (1 + fcGHz) + 0.41;
766  table3gpp->m_uLgZSD = -0.04 * log10 (1 + fcGHz) + 0.92;
767  table3gpp->m_sigLgZSD = -0.07 * log10 (1 + fcGHz) + 0.41;
768  table3gpp->m_offsetZOD = 0;
769  table3gpp->m_cDS = 11;
770  table3gpp->m_cASD = 22;
771  table3gpp->m_cASA = 22;
772  table3gpp->m_cZSA = 7;
773  table3gpp->m_uK = 0; // N/A
774  table3gpp->m_sigK = 0; // N/A
775  table3gpp->m_rTau = 2.1;
776  table3gpp->m_uXpr = 8;
777  table3gpp->m_sigXpr = 3;
778  table3gpp->m_perClusterShadowingStd = 4;
779 
780  for (uint8_t row = 0; row < 6; row++)
781  {
782  for (uint8_t column = 0; column < 6; column++)
783  {
784  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_NLOS[row][column];
785  }
786  }
787  }
788  else if (channelCondition->IsNlosv ())
789  {
790  table3gpp->m_numOfCluster = 19;
791  table3gpp->m_raysPerCluster = 20;
792  table3gpp->m_uLgDS = -0.4 * log10 (1 + fcGHz) - 7;
793  table3gpp->m_sigLgDS = 0.1;
794  table3gpp->m_uLgASD = -0.1 * log10 (1 + fcGHz) + 1.7;
795  table3gpp->m_sigLgASD = 0.1;
796  table3gpp->m_uLgASA = -0.1 * log10 (1 + fcGHz) + 1.7;
797  table3gpp->m_sigLgASA = 0.1;
798  table3gpp->m_uLgZSA = -0.04 * log10 (1 + fcGHz) + 0.92;
799  table3gpp->m_sigLgZSA = -0.07 * log10 (1 + fcGHz) + 0.41;
800  table3gpp->m_uLgZSD = -0.04 * log10 (1 + fcGHz) + 0.92;
801  table3gpp->m_sigLgZSD = -0.07 * log10 (1 + fcGHz) + 0.41;
802  table3gpp->m_offsetZOD = 0;
803  table3gpp->m_cDS = 11;
804  table3gpp->m_cASD = 22;
805  table3gpp->m_cASA = 22;
806  table3gpp->m_cZSA = 7;
807  table3gpp->m_uK = 0;
808  table3gpp->m_sigK = 4.5;
809  table3gpp->m_rTau = 2.1;
810  table3gpp->m_uXpr = 8;
811  table3gpp->m_sigXpr = 3;
812  table3gpp->m_perClusterShadowingStd = 4;
813 
814  for (uint8_t row = 0; row < 6; row++)
815  {
816  for (uint8_t column = 0; column < 6; column++)
817  {
818  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
819  }
820  }
821  }
822  else
823  {
824  NS_FATAL_ERROR ("Unknown channel condition");
825  }
826  }
827  else if (m_scenario == "V2V-Highway")
828  {
829  if (channelCondition->IsLos ())
830  {
831  table3gpp->m_numOfCluster = 12;
832  table3gpp->m_raysPerCluster = 20;
833  table3gpp->m_uLgDS = -8.3;
834  table3gpp->m_sigLgDS = 0.2;
835  table3gpp->m_uLgASD = 1.4;
836  table3gpp->m_sigLgASD = 0.1;
837  table3gpp->m_uLgASA = 1.4;
838  table3gpp->m_sigLgASA = 0.1;
839  table3gpp->m_uLgZSA = -0.1 * log10 (1 + fcGHz) + 0.73;
840  table3gpp->m_sigLgZSA = -0.04 * log10 (1 + fcGHz) + 0.34;
841  table3gpp->m_uLgZSD = -0.1 * log10 (1 + fcGHz) + 0.73;
842  table3gpp->m_sigLgZSD = -0.04 * log10 (1 + fcGHz) + 0.34;
843  table3gpp->m_offsetZOD = 0;
844  table3gpp->m_cDS = 5;
845  table3gpp->m_cASD = 17;
846  table3gpp->m_cASA = 17;
847  table3gpp->m_cZSA = 7;
848  table3gpp->m_uK = 9;
849  table3gpp->m_sigK = 3.5;
850  table3gpp->m_rTau = 3;
851  table3gpp->m_uXpr = 9;
852  table3gpp->m_sigXpr = 3;
853  table3gpp->m_perClusterShadowingStd = 4;
854 
855  for (uint8_t row = 0; row < 7; row++)
856  {
857  for (uint8_t column = 0; column < 7; column++)
858  {
859  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
860  }
861  }
862  }
863  else if (channelCondition->IsNlosv ())
864  {
865  table3gpp->m_numOfCluster = 19;
866  table3gpp->m_raysPerCluster = 20;
867  table3gpp->m_uLgDS = -8.3;
868  table3gpp->m_sigLgDS = 0.3;
869  table3gpp->m_uLgASD = 1.5;
870  table3gpp->m_sigLgASD = 0.1;
871  table3gpp->m_uLgASA = 1.5;
872  table3gpp->m_sigLgASA = 0.1;
873  table3gpp->m_uLgZSA = -0.04 * log10 (1 + fcGHz) + 0.92;
874  table3gpp->m_sigLgZSA = -0.07 * log10 (1 + fcGHz) + 0.41;
875  table3gpp->m_uLgZSD = -0.04 * log10 (1 + fcGHz) + 0.92;
876  table3gpp->m_sigLgZSD = -0.07 * log10 (1 + fcGHz) + 0.41;
877  table3gpp->m_offsetZOD = 0;
878  table3gpp->m_cDS = 11;
879  table3gpp->m_cASD = 22;
880  table3gpp->m_cASA = 22;
881  table3gpp->m_cZSA = 7;
882  table3gpp->m_uK = 0;
883  table3gpp->m_sigK = 4.5;
884  table3gpp->m_rTau = 2.1;
885  table3gpp->m_uXpr = 8.0;
886  table3gpp->m_sigXpr = 3;
887  table3gpp->m_perClusterShadowingStd = 4;
888 
889  for (uint8_t row = 0; row < 6; row++)
890  {
891  for (uint8_t column = 0; column < 6; column++)
892  {
893  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
894  }
895  }
896  }
897  else if (channelCondition->IsNlos ())
898  {
899  NS_LOG_WARN ("The fast fading parameters for the NLOS condition in the Highway scenario are not defined in TR 37.885, use the ones defined in TDoc R1-1803671 instead");
900 
901  table3gpp->m_numOfCluster = 19;
902  table3gpp->m_raysPerCluster = 20;
903  table3gpp->m_uLgDS = -0.3 * log10 (1 + fcGHz) - 7;
904  table3gpp->m_sigLgDS = 0.28;
905  table3gpp->m_uLgASD = -0.08 * log10 (1 + fcGHz) + 1.81;
906  table3gpp->m_sigLgASD = 0.05 * log10 (1 + fcGHz) + 0.3;
907  table3gpp->m_uLgASA = -0.08 * log10 (1 + fcGHz) + 1.81;
908  table3gpp->m_sigLgASA = 0.05 * log10 (1 + fcGHz) + 0.3;
909  table3gpp->m_uLgZSA = -0.04 * log10 (1 + fcGHz) + 0.92,
910  table3gpp->m_sigLgZSA = -0.07 * log10 (1 + fcGHz) + 0.41;
911  table3gpp->m_uLgZSD = -0.04 * log10 (1 + fcGHz) + 0.92;
912  table3gpp->m_sigLgZSD = -0.07 * log10 (1 + fcGHz) + 0.41;
913  table3gpp->m_offsetZOD = 0;
914  table3gpp->m_cDS = 11;
915  table3gpp->m_cASD = 22;
916  table3gpp->m_cASA = 22;
917  table3gpp->m_cZSA = 7;
918  table3gpp->m_uK = 0; // N/A
919  table3gpp->m_sigK = 0; // N/A
920  table3gpp->m_rTau = 2.1;
921  table3gpp->m_uXpr = 8;
922  table3gpp->m_sigXpr = 3;
923  table3gpp->m_perClusterShadowingStd = 4;
924 
925  for (uint8_t row = 0; row < 6; row++)
926  {
927  for (uint8_t column = 0; column < 6; column++)
928  {
929  table3gpp->m_sqrtC[row][column] = sqrtC_UMi_NLOS[row][column];
930  }
931  }
932  }
933  else
934  {
935  NS_FATAL_ERROR ("Unknown channel condition");
936  }
937  }
938  else
939  {
940  NS_FATAL_ERROR ("unkonw scenarios");
941  }
942 
943  return table3gpp;
944 }
945 
946 bool
947 ThreeGppChannelModel::ChannelMatrixNeedsUpdate (Ptr<const ThreeGppChannelMatrix> channelMatrix, Ptr<const ChannelCondition> channelCondition) const
948 {
949  NS_LOG_FUNCTION (this);
950 
951  bool update = false;
952 
953  // if the channel condition is different the channel has to be updated
954  if (!channelMatrix->m_channelCondition->IsEqual (channelCondition))
955  {
956  NS_LOG_DEBUG ("Update the channel condition");
957  update = true;
958  }
959 
960  // if the coherence time is over the channel has to be updated
961  if (!m_updatePeriod.IsZero () && Simulator::Now () - channelMatrix->m_generatedTime > m_updatePeriod)
962  {
963  NS_LOG_DEBUG ("Generation time " << channelMatrix->m_generatedTime.As (Time::NS) << " now " << Now ().As (Time::NS));
964  update = true;
965  }
966 
967  return update;
968 }
969 
970 Ptr<const MatrixBasedChannelModel::ChannelMatrix>
971 ThreeGppChannelModel::GetChannel (Ptr<const MobilityModel> aMob,
972  Ptr<const MobilityModel> bMob,
973  Ptr<const ThreeGppAntennaArrayModel> aAntenna,
974  Ptr<const ThreeGppAntennaArrayModel> bAntenna)
975 {
976  NS_LOG_FUNCTION (this);
977 
978  // Compute the channel key. The key is reciprocal, i.e., key (a, b) = key (b, a)
979  uint32_t x1 = std::min (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
980  uint32_t x2 = std::max (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
981  uint32_t channelId = GetKey (x1, x2);
982 
983  // retrieve the channel condition
984  Ptr<const ChannelCondition> condition = m_channelConditionModel->GetChannelCondition (aMob, bMob);
985 
986  // Check if the channel is present in the map and return it, otherwise
987  // generate a new channel
988  bool update = false;
989  bool notFound = false;
990  Ptr<ThreeGppChannelMatrix> channelMatrix;
991  if (m_channelMap.find (channelId) != m_channelMap.end ())
992  {
993  // channel matrix present in the map
994  NS_LOG_DEBUG ("channel matrix present in the map");
995  channelMatrix = m_channelMap[channelId];
996 
997  // check if it has to be updated
998  update = ChannelMatrixNeedsUpdate (channelMatrix, condition);
999  }
1000  else
1001  {
1002  NS_LOG_DEBUG ("channel matrix not found");
1003  notFound = true;
1004  }
1005 
1006  // If the channel is not present in the map or if it has to be updated
1007  // generate a new realization
1008  if (notFound || update)
1009  {
1010  // channel matrix not found or has to be updated, generate a new one
1011  Angles txAngle (bMob->GetPosition (), aMob->GetPosition ());
1012  Angles rxAngle (aMob->GetPosition (), bMob->GetPosition ());
1013 
1014  double x = aMob->GetPosition ().x - bMob->GetPosition ().x;
1015  double y = aMob->GetPosition ().y - bMob->GetPosition ().y;
1016  double distance2D = sqrt (x * x + y * y);
1017 
1018  // NOTE we assume hUT = min (height(a), height(b)) and
1019  // hBS = max (height (a), height (b))
1020  double hUt = std::min (aMob->GetPosition ().z, bMob->GetPosition ().z);
1021  double hBs = std::max (aMob->GetPosition ().z, bMob->GetPosition ().z);
1022 
1023  // TODO this is not currently used, it is needed for the computation of the
1024  // additional blockage in case of spatial consistent update
1025  // I do not know who is the UT, I can use the relative distance between
1026  // tx and rx instead
1027  Vector locUt = Vector (0.0, 0.0, 0.0);
1028 
1029  channelMatrix = GetNewChannel (locUt, condition, aAntenna, bAntenna, rxAngle, txAngle, distance2D, hBs, hUt);
1030  channelMatrix->m_nodeIds = std::make_pair (aMob->GetObject<Node> ()->GetId (), bMob->GetObject<Node> ()->GetId ());
1031 
1032  // store or replace the channel matrix in the channel map
1033  m_channelMap[channelId] = channelMatrix;
1034  }
1035 
1036  return channelMatrix;
1037 }
1038 
1039 Ptr<ThreeGppChannelModel::ThreeGppChannelMatrix>
1040 ThreeGppChannelModel::GetNewChannel (Vector locUT, Ptr<const ChannelCondition> channelCondition,
1041  Ptr<const ThreeGppAntennaArrayModel> sAntenna,
1042  Ptr<const ThreeGppAntennaArrayModel> uAntenna,
1043  Angles &uAngle, Angles &sAngle,
1044  double dis2D, double hBS, double hUT) const
1045 {
1046  NS_LOG_FUNCTION (this);
1047 
1048  NS_ASSERT_MSG (m_frequency > 0.0, "Set the operating frequency first!");
1049 
1050  // get the 3GPP parameters
1051  Ptr<const ParamsTable> table3gpp = GetThreeGppTable (channelCondition, hBS, hUT, dis2D);
1052 
1053  // get the number of clusters and the number of rays per cluster
1054  uint8_t numOfCluster = table3gpp->m_numOfCluster;
1055  uint8_t raysPerCluster = table3gpp->m_raysPerCluster;
1056 
1057  // create a channel matrix instance
1058  Ptr<ThreeGppChannelMatrix> channelParams = Create<ThreeGppChannelMatrix> ();
1059  channelParams->m_channelCondition = channelCondition; // set the channel condition
1060  channelParams->m_generatedTime = Simulator::Now ();
1061 
1062  // compute the 3D distance using eq. 7.4-1
1063  double dis3D = std::sqrt (dis2D * dis2D + (hBS - hUT) * (hBS - hUT));
1064 
1065  bool los = channelCondition->IsLos ();
1066  bool o2i = channelCondition->IsO2i ();
1067 
1068  //Step 4: Generate large scale parameters. All LSPS are uncorrelated.
1069  DoubleVector LSPsIndep, LSPs;
1070  uint8_t paramNum;
1071  if (los)
1072  {
1073  paramNum = 7;
1074  }
1075  else
1076  {
1077  paramNum = 6;
1078  }
1079  //Generate paramNum independent LSPs.
1080  for (uint8_t iter = 0; iter < paramNum; iter++)
1081  {
1082  LSPsIndep.push_back (m_normalRv->GetValue ());
1083  }
1084  for (uint8_t row = 0; row < paramNum; row++)
1085  {
1086  double temp = 0;
1087  for (uint8_t column = 0; column < paramNum; column++)
1088  {
1089  temp += table3gpp->m_sqrtC[row][column] * LSPsIndep[column];
1090  }
1091  LSPs.push_back (temp);
1092  }
1093 
1094  // NOTE the shadowing is generated in the propagation loss model
1095 
1096  double DS,ASD,ASA,ZSA,ZSD,K_factor = 0;
1097  if (los)
1098  {
1099  K_factor = LSPs[1] * table3gpp->m_sigK + table3gpp->m_uK;
1100  DS = pow (10, LSPs[2] * table3gpp->m_sigLgDS + table3gpp->m_uLgDS);
1101  ASD = pow (10, LSPs[3] * table3gpp->m_sigLgASD + table3gpp->m_uLgASD);
1102  ASA = pow (10, LSPs[4] * table3gpp->m_sigLgASA + table3gpp->m_uLgASA);
1103  ZSD = pow (10, LSPs[5] * table3gpp->m_sigLgZSD + table3gpp->m_uLgZSD);
1104  ZSA = pow (10, LSPs[6] * table3gpp->m_sigLgZSA + table3gpp->m_uLgZSA);
1105  }
1106  else
1107  {
1108  DS = pow (10, LSPs[1] * table3gpp->m_sigLgDS + table3gpp->m_uLgDS);
1109  ASD = pow (10, LSPs[2] * table3gpp->m_sigLgASD + table3gpp->m_uLgASD);
1110  ASA = pow (10, LSPs[3] * table3gpp->m_sigLgASA + table3gpp->m_uLgASA);
1111  ZSD = pow (10, LSPs[4] * table3gpp->m_sigLgZSD + table3gpp->m_uLgZSD);
1112  ZSA = pow (10, LSPs[5] * table3gpp->m_sigLgZSA + table3gpp->m_uLgZSA);
1113 
1114  }
1115  ASD = std::min (ASD, 104.0);
1116  ASA = std::min (ASA, 104.0);
1117  ZSD = std::min (ZSD, 52.0);
1118  ZSA = std::min (ZSA, 52.0);
1119 
1120  channelParams->m_DS = DS;
1121  channelParams->m_K = K_factor;
1122 
1123  NS_LOG_INFO ("K-factor=" << K_factor << ",DS=" << DS << ", ASD=" << ASD << ", ASA=" << ASA << ", ZSD=" << ZSD << ", ZSA=" << ZSA);
1124 
1125  //Step 5: Generate Delays.
1126  DoubleVector clusterDelay;
1127  double minTau = 100.0;
1128  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1129  {
1130  double tau = -1 * table3gpp->m_rTau * DS * log (m_uniformRv->GetValue (0,1)); //(7.5-1)
1131  if (minTau > tau)
1132  {
1133  minTau = tau;
1134  }
1135  clusterDelay.push_back (tau);
1136  }
1137 
1138  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1139  {
1140  clusterDelay[cIndex] -= minTau;
1141  }
1142  std::sort (clusterDelay.begin (), clusterDelay.end ()); //(7.5-2)
1143 
1144  /* since the scaled Los delays are not to be used in cluster power generation,
1145  * we will generate cluster power first and resume to compute Los cluster delay later.*/
1146 
1147  //Step 6: Generate cluster powers.
1148  DoubleVector clusterPower;
1149  double powerSum = 0;
1150  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1151  {
1152  double power = exp (-1 * clusterDelay[cIndex] * (table3gpp->m_rTau - 1) / table3gpp->m_rTau / DS) *
1153  pow (10,-1 * m_normalRv->GetValue () * table3gpp->m_perClusterShadowingStd / 10); //(7.5-5)
1154  powerSum += power;
1155  clusterPower.push_back (power);
1156  }
1157  double powerMax = 0;
1158 
1159  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1160  {
1161  clusterPower[cIndex] = clusterPower[cIndex] / powerSum; //(7.5-6)
1162  }
1163 
1164  DoubleVector clusterPowerForAngles; // this power is only for equation (7.5-9) and (7.5-14), not for (7.5-22)
1165  if (los)
1166  {
1167  double K_linear = pow (10,K_factor / 10);
1168 
1169  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1170  {
1171  if (cIndex == 0)
1172  {
1173  clusterPowerForAngles.push_back (clusterPower[cIndex] / (1 + K_linear) + K_linear / (1 + K_linear)); //(7.5-8)
1174  }
1175  else
1176  {
1177  clusterPowerForAngles.push_back (clusterPower[cIndex] / (1 + K_linear)); //(7.5-8)
1178  }
1179  if (powerMax < clusterPowerForAngles[cIndex])
1180  {
1181  powerMax = clusterPowerForAngles[cIndex];
1182  }
1183  }
1184  }
1185  else
1186  {
1187  for (uint8_t cIndex = 0; cIndex < numOfCluster; cIndex++)
1188  {
1189  clusterPowerForAngles.push_back (clusterPower[cIndex]); //(7.5-6)
1190  if (powerMax < clusterPowerForAngles[cIndex])
1191  {
1192  powerMax = clusterPowerForAngles[cIndex];
1193  }
1194  }
1195  }
1196 
1197  //remove clusters with less than -25 dB power compared to the maxim cluster power;
1198  //double thresh = pow(10,-2.5);
1199  double thresh = 0.0032;
1200  for (uint8_t cIndex = numOfCluster; cIndex > 0; cIndex--)
1201  {
1202  if (clusterPowerForAngles[cIndex - 1] < thresh * powerMax )
1203  {
1204  clusterPowerForAngles.erase (clusterPowerForAngles.begin () + cIndex - 1);
1205  clusterPower.erase (clusterPower.begin () + cIndex - 1);
1206  clusterDelay.erase (clusterDelay.begin () + cIndex - 1);
1207  }
1208  }
1209 
1210  NS_ASSERT (clusterPower.size () < UINT8_MAX);
1211  uint8_t numReducedCluster = clusterPower.size ();
1212 
1213  channelParams->m_numCluster = numReducedCluster;
1214  // Resume step 5 to compute the delay for LoS condition.
1215  if (los)
1216  {
1217  double C_tau = 0.7705 - 0.0433 * K_factor + 2e-4 * pow (K_factor,2) + 17e-6 * pow (K_factor,3); //(7.5-3)
1218  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1219  {
1220  clusterDelay[cIndex] = clusterDelay[cIndex] / C_tau; //(7.5-4)
1221  }
1222  }
1223 
1224  //Step 7: Generate arrival and departure angles for both azimuth and elevation.
1225 
1226  double C_NLOS, C_phi;
1227  //According to table 7.5-6, only cluster number equals to 8, 10, 11, 12, 19 and 20 is valid.
1228  //Not sure why the other cases are in Table 7.5-2.
1229  switch (numOfCluster) // Table 7.5-2
1230  {
1231  case 4:
1232  C_NLOS = 0.779;
1233  break;
1234  case 5:
1235  C_NLOS = 0.860;
1236  break;
1237  case 8:
1238  C_NLOS = 1.018;
1239  break;
1240  case 10:
1241  C_NLOS = 1.090;
1242  break;
1243  case 11:
1244  C_NLOS = 1.123;
1245  break;
1246  case 12:
1247  C_NLOS = 1.146;
1248  break;
1249  case 14:
1250  C_NLOS = 1.190;
1251  break;
1252  case 15:
1253  C_NLOS = 1.221;
1254  break;
1255  case 16:
1256  C_NLOS = 1.226;
1257  break;
1258  case 19:
1259  C_NLOS = 1.273;
1260  break;
1261  case 20:
1262  C_NLOS = 1.289;
1263  break;
1264  default:
1265  NS_FATAL_ERROR ("Invalid cluster number");
1266  }
1267 
1268  if (los)
1269  {
1270  C_phi = C_NLOS * (1.1035 - 0.028 * K_factor - 2e-3 * pow (K_factor,2) + 1e-4 * pow (K_factor,3)); //(7.5-10))
1271  }
1272  else
1273  {
1274  C_phi = C_NLOS; //(7.5-10)
1275  }
1276 
1277  double C_theta;
1278  switch (numOfCluster) //Table 7.5-4
1279  {
1280  case 8:
1281  C_NLOS = 0.889;
1282  break;
1283  case 10:
1284  C_NLOS = 0.957;
1285  break;
1286  case 11:
1287  C_NLOS = 1.031;
1288  break;
1289  case 12:
1290  C_NLOS = 1.104;
1291  break;
1292  case 15:
1293  C_NLOS = 1.1088;
1294  break;
1295  case 19:
1296  C_NLOS = 1.184;
1297  break;
1298  case 20:
1299  C_NLOS = 1.178;
1300  break;
1301  default:
1302  NS_FATAL_ERROR ("Invalid cluster number");
1303  }
1304 
1305  if (los)
1306  {
1307  C_theta = C_NLOS * (1.3086 + 0.0339 * K_factor - 0.0077 * pow (K_factor,2) + 2e-4 * pow (K_factor,3)); //(7.5-15)
1308  }
1309  else
1310  {
1311  C_theta = C_NLOS;
1312  }
1313 
1314 
1315  DoubleVector clusterAoa, clusterAod, clusterZoa, clusterZod;
1316  double angle;
1317  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1318  {
1319  angle = 2 * ASA * sqrt (-1 * log (clusterPowerForAngles[cIndex] / powerMax)) / 1.4 / C_phi; //(7.5-9)
1320  clusterAoa.push_back (angle);
1321  angle = 2 * ASD * sqrt (-1 * log (clusterPowerForAngles[cIndex] / powerMax)) / 1.4 / C_phi; //(7.5-9)
1322  clusterAod.push_back (angle);
1323  angle = -1 * ZSA * log (clusterPowerForAngles[cIndex] / powerMax) / C_theta; //(7.5-14)
1324  clusterZoa.push_back (angle);
1325  angle = -1 * ZSD * log (clusterPowerForAngles[cIndex] / powerMax) / C_theta;
1326  clusterZod.push_back (angle);
1327  }
1328 
1329  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1330  {
1331  int Xn = 1;
1332  if (m_uniformRv->GetValue (0,1) < 0.5)
1333  {
1334  Xn = -1;
1335  }
1336  clusterAoa[cIndex] = clusterAoa[cIndex] * Xn + (m_normalRv->GetValue () * ASA / 7) + uAngle.phi * 180 / M_PI; //(7.5-11)
1337  clusterAod[cIndex] = clusterAod[cIndex] * Xn + (m_normalRv->GetValue () * ASD / 7) + sAngle.phi * 180 / M_PI;
1338  if (o2i)
1339  {
1340  clusterZoa[cIndex] = clusterZoa[cIndex] * Xn + (m_normalRv->GetValue () * ZSA / 7) + 90; //(7.5-16)
1341  }
1342  else
1343  {
1344  clusterZoa[cIndex] = clusterZoa[cIndex] * Xn + (m_normalRv->GetValue () * ZSA / 7) + uAngle.theta * 180 / M_PI; //(7.5-16)
1345  }
1346  clusterZod[cIndex] = clusterZod[cIndex] * Xn + (m_normalRv->GetValue () * ZSD / 7) + sAngle.theta * 180 / M_PI + table3gpp->m_offsetZOD; //(7.5-19)
1347 
1348  }
1349 
1350  if (los)
1351  {
1352  //The 7.5-12 can be rewrite as Theta_n,ZOA = Theta_n,ZOA - (Theta_1,ZOA - Theta_LOS,ZOA) = Theta_n,ZOA - diffZOA,
1353  //Similar as AOD, ZSA and ZSD.
1354  double diffAoa = clusterAoa[0] - uAngle.phi * 180 / M_PI;
1355  double diffAod = clusterAod[0] - sAngle.phi * 180 / M_PI;
1356  double diffZsa = clusterZoa[0] - uAngle.theta * 180 / M_PI;
1357  double diffZsd = clusterZod[0] - sAngle.theta * 180 / M_PI;
1358 
1359  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1360  {
1361  clusterAoa[cIndex] -= diffAoa; //(7.5-12)
1362  clusterAod[cIndex] -= diffAod;
1363  clusterZoa[cIndex] -= diffZsa; //(7.5-17)
1364  clusterZod[cIndex] -= diffZsd;
1365 
1366  }
1367  }
1368 
1369  double rayAoa_radian[numReducedCluster][raysPerCluster]; //rayAoa_radian[n][m], where n is cluster index, m is ray index
1370  double rayAod_radian[numReducedCluster][raysPerCluster]; //rayAod_radian[n][m], where n is cluster index, m is ray index
1371  double rayZoa_radian[numReducedCluster][raysPerCluster]; //rayZoa_radian[n][m], where n is cluster index, m is ray index
1372  double rayZod_radian[numReducedCluster][raysPerCluster]; //rayZod_radian[n][m], where n is cluster index, m is ray index
1373 
1374  for (uint8_t nInd = 0; nInd < numReducedCluster; nInd++)
1375  {
1376  for (uint8_t mInd = 0; mInd < raysPerCluster; mInd++)
1377  {
1378  double tempAoa = clusterAoa[nInd] + table3gpp->m_cASA * offSetAlpha[mInd]; //(7.5-13)
1379  while (tempAoa > 360)
1380  {
1381  tempAoa -= 360;
1382  }
1383 
1384  while (tempAoa < 0)
1385  {
1386  tempAoa += 360;
1387 
1388  }
1389  NS_ASSERT_MSG (tempAoa >= 0 && tempAoa <= 360, "the AOA should be the range of [0,360]");
1390  rayAoa_radian[nInd][mInd] = tempAoa * M_PI / 180;
1391 
1392  double tempAod = clusterAod[nInd] + table3gpp->m_cASD * offSetAlpha[mInd];
1393  while (tempAod > 360)
1394  {
1395  tempAod -= 360;
1396  }
1397 
1398  while (tempAod < 0)
1399  {
1400  tempAod += 360;
1401  }
1402  NS_ASSERT_MSG (tempAod >= 0 && tempAod <= 360, "the AOD should be the range of [0,360]");
1403  rayAod_radian[nInd][mInd] = tempAod * M_PI / 180;
1404 
1405  double tempZoa = clusterZoa[nInd] + table3gpp->m_cZSA * offSetAlpha[mInd]; //(7.5-18)
1406 
1407  while (tempZoa > 360)
1408  {
1409  tempZoa -= 360;
1410  }
1411 
1412  while (tempZoa < 0)
1413  {
1414  tempZoa += 360;
1415  }
1416 
1417  if (tempZoa > 180)
1418  {
1419  tempZoa = 360 - tempZoa;
1420  }
1421 
1422  NS_ASSERT_MSG (tempZoa >= 0&&tempZoa <= 180, "the ZOA should be the range of [0,180]");
1423  rayZoa_radian[nInd][mInd] = tempZoa * M_PI / 180;
1424 
1425  double tempZod = clusterZod[nInd] + 0.375 * pow (10,table3gpp->m_uLgZSD) * offSetAlpha[mInd]; //(7.5-20)
1426 
1427  while (tempZod > 360)
1428  {
1429  tempZod -= 360;
1430  }
1431 
1432  while (tempZod < 0)
1433  {
1434  tempZod += 360;
1435  }
1436  if (tempZod > 180)
1437  {
1438  tempZod = 360 - tempZod;
1439  }
1440  NS_ASSERT_MSG (tempZod >= 0&&tempZod <= 180, "the ZOD should be the range of [0,180]");
1441  rayZod_radian[nInd][mInd] = tempZod * M_PI / 180;
1442  }
1443  }
1444  DoubleVector angle_degree;
1445  double sizeTemp = clusterZoa.size ();
1446  for (uint8_t ind = 0; ind < 4; ind++)
1447  {
1448  switch (ind)
1449  {
1450  case 0:
1451  angle_degree = clusterAoa;
1452  break;
1453  case 1:
1454  angle_degree = clusterZoa;
1455  break;
1456  case 2:
1457  angle_degree = clusterAod;
1458  break;
1459  case 3:
1460  angle_degree = clusterZod;
1461  break;
1462  default:
1463  NS_FATAL_ERROR ("Programming Error");
1464  }
1465 
1466  for (uint8_t nIndex = 0; nIndex < sizeTemp; nIndex++)
1467  {
1468  while (angle_degree[nIndex] > 360)
1469  {
1470  angle_degree[nIndex] -= 360;
1471  }
1472 
1473  while (angle_degree[nIndex] < 0)
1474  {
1475  angle_degree[nIndex] += 360;
1476  }
1477 
1478  if (ind == 1 || ind == 3)
1479  {
1480  if (angle_degree[nIndex] > 180)
1481  {
1482  angle_degree[nIndex] = 360 - angle_degree[nIndex];
1483  }
1484  }
1485  }
1486  switch (ind)
1487  {
1488  case 0:
1489  clusterAoa = angle_degree;
1490  break;
1491  case 1:
1492  clusterZoa = angle_degree;
1493  break;
1494  case 2:
1495  clusterAod = angle_degree;
1496  break;
1497  case 3:
1498  clusterZod = angle_degree;
1499  break;
1500  default:
1501  NS_FATAL_ERROR ("Programming Error");
1502  }
1503  }
1504 
1505  DoubleVector attenuation_dB;
1506  if (m_blockage)
1507  {
1508  attenuation_dB = CalcAttenuationOfBlockage (channelParams, clusterAoa, clusterZoa);
1509  for (uint8_t cInd = 0; cInd < numReducedCluster; cInd++)
1510  {
1511  clusterPower[cInd] = clusterPower[cInd] / pow (10,attenuation_dB[cInd] / 10);
1512  }
1513  }
1514  else
1515  {
1516  attenuation_dB.push_back (0);
1517  }
1518 
1519  //Step 8: Coupling of rays within a cluster for both azimuth and elevation
1520  //shuffle all the arrays to perform random coupling
1521  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1522  {
1523  Shuffle (&rayAod_radian[cIndex][0], &rayAod_radian[cIndex][raysPerCluster]);
1524  Shuffle (&rayAoa_radian[cIndex][0], &rayAoa_radian[cIndex][raysPerCluster]);
1525  Shuffle (&rayZod_radian[cIndex][0], &rayZod_radian[cIndex][raysPerCluster]);
1526  Shuffle (&rayZoa_radian[cIndex][0], &rayZoa_radian[cIndex][raysPerCluster]);
1527  }
1528 
1529  //Step 9: Generate the cross polarization power ratios
1530  //Step 10: Draw initial phases
1531  Double2DVector crossPolarizationPowerRatios; // vector containing the cross polarization power ratios, as defined by 7.5-21
1532  Double3DVector clusterPhase; //rayAoa_radian[n][m], where n is cluster index, m is ray index
1533  for (uint8_t nInd = 0; nInd < numReducedCluster; nInd++)
1534  {
1535  DoubleVector temp; // used to store the XPR values
1536  Double2DVector temp2; // used to store the PHI values for all the possible combination of polarization
1537  for (uint8_t mInd = 0; mInd < raysPerCluster; mInd++)
1538  {
1539  double uXprLinear = pow (10, table3gpp->m_uXpr / 10); // convert to linear
1540  double sigXprLinear = pow (10, table3gpp->m_sigXpr / 10); // convert to linear
1541 
1542  temp.push_back (std::pow (10, (m_normalRv->GetValue () * sigXprLinear + uXprLinear) / 10));
1543  DoubleVector temp3; // used to store the PHI valuse
1544  for (uint8_t pInd = 0; pInd < 4; pInd++)
1545  {
1546  temp3.push_back (m_uniformRv->GetValue (-1 * M_PI, M_PI));
1547  }
1548  temp2.push_back (temp3);
1549  }
1550  crossPolarizationPowerRatios.push_back (temp);
1551  clusterPhase.push_back (temp2);
1552  }
1553  channelParams->m_clusterPhase = clusterPhase;
1554 
1555  //Step 11: Generate channel coefficients for each cluster n and each receiver
1556  // and transmitter element pair u,s.
1557 
1558  Complex3DVector H_NLOS; // channel coefficients H_NLOS [u][s][n],
1559  // where u and s are receive and transmit antenna element, n is cluster index.
1560  uint64_t uSize = uAntenna->GetNumberOfElements ();
1561  uint64_t sSize = sAntenna->GetNumberOfElements ();
1562 
1563  uint8_t cluster1st = 0, cluster2nd = 0; // first and second strongest cluster;
1564  double maxPower = 0;
1565  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1566  {
1567  if (maxPower < clusterPower[cIndex])
1568  {
1569  maxPower = clusterPower[cIndex];
1570  cluster1st = cIndex;
1571  }
1572  }
1573  maxPower = 0;
1574  for (uint8_t cIndex = 0; cIndex < numReducedCluster; cIndex++)
1575  {
1576  if (maxPower < clusterPower[cIndex] && cluster1st != cIndex)
1577  {
1578  maxPower = clusterPower[cIndex];
1579  cluster2nd = cIndex;
1580  }
1581  }
1582 
1583  NS_LOG_INFO ("1st strongest cluster:" << (int)cluster1st << ", 2nd strongest cluster:" << (int)cluster2nd);
1584 
1585  Complex3DVector H_usn; //channel coffecient H_usn[u][s][n];
1586  // NOTE Since each of the strongest 2 clusters are divided into 3 sub-clusters,
1587  // the total cluster will be numReducedCLuster + 4.
1588 
1589  H_usn.resize (uSize);
1590  for (uint64_t uIndex = 0; uIndex < uSize; uIndex++)
1591  {
1592  H_usn[uIndex].resize (sSize);
1593  for (uint64_t sIndex = 0; sIndex < sSize; sIndex++)
1594  {
1595  H_usn[uIndex][sIndex].resize (numReducedCluster);
1596  }
1597  }
1598 
1599  // The following for loops computes the channel coefficients
1600  for (uint64_t uIndex = 0; uIndex < uSize; uIndex++)
1601  {
1602  Vector uLoc = uAntenna->GetElementLocation (uIndex);
1603 
1604  for (uint64_t sIndex = 0; sIndex < sSize; sIndex++)
1605  {
1606 
1607  Vector sLoc = sAntenna->GetElementLocation (sIndex);
1608 
1609  for (uint8_t nIndex = 0; nIndex < numReducedCluster; nIndex++)
1610  {
1611  //Compute the N-2 weakest cluster, only vertical polarization. (7.5-22)
1612  if (nIndex != cluster1st && nIndex != cluster2nd)
1613  {
1614  std::complex<double> rays (0,0);
1615  for (uint8_t mIndex = 0; mIndex < raysPerCluster; mIndex++)
1616  {
1617  DoubleVector initialPhase = clusterPhase[nIndex][mIndex];
1618  double k = crossPolarizationPowerRatios[nIndex][mIndex];
1619  //lambda_0 is accounted in the antenna spacing uLoc and sLoc.
1620  double rxPhaseDiff = 2 * M_PI * (sin (rayZoa_radian[nIndex][mIndex]) * cos (rayAoa_radian[nIndex][mIndex]) * uLoc.x
1621  + sin (rayZoa_radian[nIndex][mIndex]) * sin (rayAoa_radian[nIndex][mIndex]) * uLoc.y
1622  + cos (rayZoa_radian[nIndex][mIndex]) * uLoc.z);
1623 
1624  double txPhaseDiff = 2 * M_PI * (sin (rayZod_radian[nIndex][mIndex]) * cos (rayAod_radian[nIndex][mIndex]) * sLoc.x
1625  + sin (rayZod_radian[nIndex][mIndex]) * sin (rayAod_radian[nIndex][mIndex]) * sLoc.y
1626  + cos (rayZod_radian[nIndex][mIndex]) * sLoc.z);
1627  // NOTE Doppler is computed in the CalcBeamformingGain function and is simplified to only account for the center anngle of each cluster.
1628 
1629  double rxFieldPatternPhi, rxFieldPatternTheta, txFieldPatternPhi, txFieldPatternTheta;
1630  std::tie (rxFieldPatternPhi, rxFieldPatternTheta) = uAntenna->GetElementFieldPattern (Angles (rayAoa_radian[nIndex][mIndex], rayZoa_radian[nIndex][mIndex]));
1631  std::tie (txFieldPatternPhi, txFieldPatternTheta) = sAntenna->GetElementFieldPattern (Angles (rayAod_radian[nIndex][mIndex], rayZod_radian[nIndex][mIndex]));
1632 
1633  rays += (exp (std::complex<double> (0, initialPhase[0])) * rxFieldPatternTheta * txFieldPatternTheta +
1634  +exp (std::complex<double> (0, initialPhase[1])) * std::sqrt (1 / k) * rxFieldPatternTheta * txFieldPatternPhi +
1635  +exp (std::complex<double> (0, initialPhase[2])) * std::sqrt (1 / k) * rxFieldPatternPhi * txFieldPatternTheta +
1636  +exp (std::complex<double> (0, initialPhase[3])) * rxFieldPatternPhi * txFieldPatternPhi)
1637  * exp (std::complex<double> (0, rxPhaseDiff))
1638  * exp (std::complex<double> (0, txPhaseDiff));
1639  }
1640  rays *= sqrt (clusterPower[nIndex] / raysPerCluster);
1641  H_usn[uIndex][sIndex][nIndex] = rays;
1642  }
1643  else //(7.5-28)
1644  {
1645  std::complex<double> raysSub1 (0,0);
1646  std::complex<double> raysSub2 (0,0);
1647  std::complex<double> raysSub3 (0,0);
1648 
1649  for (uint8_t mIndex = 0; mIndex < raysPerCluster; mIndex++)
1650  {
1651  double k = crossPolarizationPowerRatios[nIndex][mIndex];
1652 
1653  //ZML:Just remind me that the angle offsets for the 3 subclusters were not generated correctly.
1654 
1655  DoubleVector initialPhase = clusterPhase[nIndex][mIndex];
1656  double rxPhaseDiff = 2 * M_PI * (sin (rayZoa_radian[nIndex][mIndex]) * cos (rayAoa_radian[nIndex][mIndex]) * uLoc.x
1657  + sin (rayZoa_radian[nIndex][mIndex]) * sin (rayAoa_radian[nIndex][mIndex]) * uLoc.y
1658  + cos (rayZoa_radian[nIndex][mIndex]) * uLoc.z);
1659  double txPhaseDiff = 2 * M_PI * (sin (rayZod_radian[nIndex][mIndex]) * cos (rayAod_radian[nIndex][mIndex]) * sLoc.x
1660  + sin (rayZod_radian[nIndex][mIndex]) * sin (rayAod_radian[nIndex][mIndex]) * sLoc.y
1661  + cos (rayZod_radian[nIndex][mIndex]) * sLoc.z);
1662 
1663  double rxFieldPatternPhi, rxFieldPatternTheta, txFieldPatternPhi, txFieldPatternTheta;
1664  std::tie (rxFieldPatternPhi, rxFieldPatternTheta) = uAntenna->GetElementFieldPattern (Angles (rayAoa_radian[nIndex][mIndex], rayZoa_radian[nIndex][mIndex]));
1665  std::tie (txFieldPatternPhi, txFieldPatternTheta) = sAntenna->GetElementFieldPattern (Angles (rayAod_radian[nIndex][mIndex], rayZod_radian[nIndex][mIndex]));
1666 
1667  switch (mIndex)
1668  {
1669  case 9:
1670  case 10:
1671  case 11:
1672  case 12:
1673  case 17:
1674  case 18:
1675  raysSub2 += (exp (std::complex<double> (0, initialPhase[0])) * rxFieldPatternTheta * txFieldPatternTheta +
1676  +exp (std::complex<double> (0, initialPhase[1])) * std::sqrt (1 / k) * rxFieldPatternTheta * txFieldPatternPhi +
1677  +exp (std::complex<double> (0, initialPhase[2])) * std::sqrt (1 / k) * rxFieldPatternPhi * txFieldPatternTheta +
1678  +exp (std::complex<double> (0, initialPhase[3])) * rxFieldPatternPhi * txFieldPatternPhi)
1679  * exp (std::complex<double> (0, rxPhaseDiff))
1680  * exp (std::complex<double> (0, txPhaseDiff));
1681  break;
1682  case 13:
1683  case 14:
1684  case 15:
1685  case 16:
1686  raysSub3 += (exp (std::complex<double> (0, initialPhase[0])) * rxFieldPatternTheta * txFieldPatternTheta +
1687  +exp (std::complex<double> (0, initialPhase[1])) * std::sqrt (1 / k) * rxFieldPatternTheta * txFieldPatternPhi +
1688  +exp (std::complex<double> (0, initialPhase[2])) * std::sqrt (1 / k) * rxFieldPatternPhi * txFieldPatternTheta +
1689  +exp (std::complex<double> (0, initialPhase[3])) * rxFieldPatternPhi * txFieldPatternPhi)
1690  * exp (std::complex<double> (0, rxPhaseDiff))
1691  * exp (std::complex<double> (0, txPhaseDiff));
1692  break;
1693  default: //case 1,2,3,4,5,6,7,8,19,20
1694  raysSub1 += (exp (std::complex<double> (0, initialPhase[0])) * rxFieldPatternTheta * txFieldPatternTheta +
1695  +exp (std::complex<double> (0, initialPhase[1])) * std::sqrt (1 / k) * rxFieldPatternTheta * txFieldPatternPhi +
1696  +exp (std::complex<double> (0, initialPhase[2])) * std::sqrt (1 / k) * rxFieldPatternPhi * txFieldPatternTheta +
1697  +exp (std::complex<double> (0, initialPhase[3])) * rxFieldPatternPhi * txFieldPatternPhi)
1698  * exp (std::complex<double> (0, rxPhaseDiff))
1699  * exp (std::complex<double> (0, txPhaseDiff));
1700  break;
1701  }
1702  }
1703  raysSub1 *= sqrt (clusterPower[nIndex] / raysPerCluster);
1704  raysSub2 *= sqrt (clusterPower[nIndex] / raysPerCluster);
1705  raysSub3 *= sqrt (clusterPower[nIndex] / raysPerCluster);
1706  H_usn[uIndex][sIndex][nIndex] = raysSub1;
1707  H_usn[uIndex][sIndex].push_back (raysSub2);
1708  H_usn[uIndex][sIndex].push_back (raysSub3);
1709 
1710  }
1711  }
1712  if (los) //(7.5-29) && (7.5-30)
1713  {
1714  std::complex<double> ray (0,0);
1715  double rxPhaseDiff = 2 * M_PI * (sin (uAngle.theta) * cos (uAngle.phi) * uLoc.x
1716  + sin (uAngle.theta) * sin (uAngle.phi) * uLoc.y
1717  + cos (uAngle.theta) * uLoc.z);
1718  double txPhaseDiff = 2 * M_PI * (sin (sAngle.theta) * cos (sAngle.phi) * sLoc.x
1719  + sin (sAngle.theta) * sin (sAngle.phi) * sLoc.y
1720  + cos (sAngle.theta) * sLoc.z);
1721 
1722  double rxFieldPatternPhi, rxFieldPatternTheta, txFieldPatternPhi, txFieldPatternTheta;
1723  std::tie (rxFieldPatternPhi, rxFieldPatternTheta) = uAntenna->GetElementFieldPattern (Angles (uAngle.phi, uAngle.theta));
1724  std::tie (txFieldPatternPhi, txFieldPatternTheta) = sAntenna->GetElementFieldPattern (Angles (sAngle.phi, sAngle.theta));
1725 
1726  double lambda = 3e8 / m_frequency; // the wavelength of the carrier frequency
1727 
1728  ray = (rxFieldPatternTheta * txFieldPatternTheta - rxFieldPatternPhi * txFieldPatternPhi)
1729  * exp (std::complex<double> (0, -2 * M_PI * dis3D / lambda))
1730  * exp (std::complex<double> (0, rxPhaseDiff))
1731  * exp (std::complex<double> (0, txPhaseDiff));
1732 
1733  double K_linear = pow (10,K_factor / 10);
1734  // the LOS path should be attenuated if blockage is enabled.
1735  H_usn[uIndex][sIndex][0] = sqrt (1 / (K_linear + 1)) * H_usn[uIndex][sIndex][0] + sqrt (K_linear / (1 + K_linear)) * ray / pow (10,attenuation_dB[0] / 10); //(7.5-30) for tau = tau1
1736  double tempSize = H_usn[uIndex][sIndex].size ();
1737  for (uint8_t nIndex = 1; nIndex < tempSize; nIndex++)
1738  {
1739  H_usn[uIndex][sIndex][nIndex] *= sqrt (1 / (K_linear + 1)); //(7.5-30) for tau = tau2...taunN
1740  }
1741 
1742  }
1743  }
1744  }
1745 
1746  // store the delays and the angles for the subclusters
1747  if (cluster1st == cluster2nd)
1748  {
1749  clusterDelay.push_back (clusterDelay[cluster1st] + 1.28 * table3gpp->m_cDS);
1750  clusterDelay.push_back (clusterDelay[cluster1st] + 2.56 * table3gpp->m_cDS);
1751 
1752  clusterAoa.push_back (clusterAoa[cluster1st]);
1753  clusterAoa.push_back (clusterAoa[cluster1st]);
1754 
1755  clusterZoa.push_back (clusterZoa[cluster1st]);
1756  clusterZoa.push_back (clusterZoa[cluster1st]);
1757 
1758  clusterAod.push_back (clusterAod[cluster1st]);
1759  clusterAod.push_back (clusterAod[cluster1st]);
1760 
1761  clusterZod.push_back (clusterZod[cluster1st]);
1762  clusterZod.push_back (clusterZod[cluster1st]);
1763  }
1764  else
1765  {
1766  double min, max;
1767  if (cluster1st < cluster2nd)
1768  {
1769  min = cluster1st;
1770  max = cluster2nd;
1771  }
1772  else
1773  {
1774  min = cluster2nd;
1775  max = cluster1st;
1776  }
1777  clusterDelay.push_back (clusterDelay[min] + 1.28 * table3gpp->m_cDS);
1778  clusterDelay.push_back (clusterDelay[min] + 2.56 * table3gpp->m_cDS);
1779  clusterDelay.push_back (clusterDelay[max] + 1.28 * table3gpp->m_cDS);
1780  clusterDelay.push_back (clusterDelay[max] + 2.56 * table3gpp->m_cDS);
1781 
1782  clusterAoa.push_back (clusterAoa[min]);
1783  clusterAoa.push_back (clusterAoa[min]);
1784  clusterAoa.push_back (clusterAoa[max]);
1785  clusterAoa.push_back (clusterAoa[max]);
1786 
1787  clusterZoa.push_back (clusterZoa[min]);
1788  clusterZoa.push_back (clusterZoa[min]);
1789  clusterZoa.push_back (clusterZoa[max]);
1790  clusterZoa.push_back (clusterZoa[max]);
1791 
1792  clusterAod.push_back (clusterAod[min]);
1793  clusterAod.push_back (clusterAod[min]);
1794  clusterAod.push_back (clusterAod[max]);
1795  clusterAod.push_back (clusterAod[max]);
1796 
1797  clusterZod.push_back (clusterZod[min]);
1798  clusterZod.push_back (clusterZod[min]);
1799  clusterZod.push_back (clusterZod[max]);
1800  clusterZod.push_back (clusterZod[max]);
1801 
1802 
1803  }
1804 
1805  NS_LOG_INFO ("size of coefficient matrix =[" << H_usn.size () << "][" << H_usn[0].size () << "][" << H_usn[0][0].size () << "]");
1806 
1807  channelParams->m_channel = H_usn;
1808  channelParams->m_delay = clusterDelay;
1809 
1810  channelParams->m_angle.clear ();
1811  channelParams->m_angle.push_back (clusterAoa);
1812  channelParams->m_angle.push_back (clusterZoa);
1813  channelParams->m_angle.push_back (clusterAod);
1814  channelParams->m_angle.push_back (clusterZod);
1815 
1816  return channelParams;
1817 }
1818 
1819 MatrixBasedChannelModel::DoubleVector
1820 ThreeGppChannelModel::CalcAttenuationOfBlockage (Ptr<ThreeGppChannelModel::ThreeGppChannelMatrix> params,
1821  const DoubleVector &clusterAOA,
1822  const DoubleVector &clusterZOA) const
1823 {
1824  NS_LOG_FUNCTION (this);
1825 
1826  DoubleVector powerAttenuation;
1827  uint8_t clusterNum = clusterAOA.size ();
1828  for (uint8_t cInd = 0; cInd < clusterNum; cInd++)
1829  {
1830  powerAttenuation.push_back (0); //Initial power attenuation for all clusters to be 0 dB;
1831  }
1832  //step a: the number of non-self blocking blockers is stored in m_numNonSelfBlocking.
1833 
1834  //step b:Generate the size and location of each blocker
1835  //generate self blocking (i.e., for blockage from the human body)
1836  double phi_sb, x_sb, theta_sb, y_sb;
1837  //table 7.6.4.1-1 Self-blocking region parameters.
1838  if (m_portraitMode)
1839  {
1840  phi_sb = 260;
1841  x_sb = 120;
1842  theta_sb = 100;
1843  y_sb = 80;
1844  }
1845  else // landscape mode
1846  {
1847  phi_sb = 40;
1848  x_sb = 160;
1849  theta_sb = 110;
1850  y_sb = 75;
1851  }
1852 
1853  //generate or update non-self blocking
1854  if (params->m_nonSelfBlocking.size () == 0) //generate new blocking regions
1855  {
1856  for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
1857  {
1858  //draw value from table 7.6.4.1-2 Blocking region parameters
1859  DoubleVector table;
1860  table.push_back (m_normalRv->GetValue ()); //phi_k: store the normal RV that will be mapped to uniform (0,360) later.
1861  if (m_scenario == "InH-OfficeMixed" || m_scenario == "InH-OfficeOpen")
1862  {
1863  table.push_back (m_uniformRv->GetValue (15, 45)); //x_k
1864  table.push_back (90); //Theta_k
1865  table.push_back (m_uniformRv->GetValue (5, 15)); //y_k
1866  table.push_back (2); //r
1867  }
1868  else
1869  {
1870  table.push_back (m_uniformRv->GetValue (5, 15)); //x_k
1871  table.push_back (90); //Theta_k
1872  table.push_back (5); //y_k
1873  table.push_back (10); //r
1874  }
1875  params->m_nonSelfBlocking.push_back (table);
1876  }
1877  }
1878  else
1879  {
1880  double deltaX = sqrt (pow (params->m_preLocUT.x - params->m_locUT.x, 2) + pow (params->m_preLocUT.y - params->m_locUT.y, 2));
1881  //if deltaX and speed are both 0, the autocorrelation is 1, skip updating
1882  if (deltaX > 1e-6 || m_blockerSpeed > 1e-6)
1883  {
1884  double corrDis;
1885  //draw value from table 7.6.4.1-4: Spatial correlation distance for different m_scenarios.
1886  if (m_scenario == "InH-OfficeMixed" || m_scenario == "InH-OfficeOpen")
1887  {
1888  //InH, correlation distance = 5;
1889  corrDis = 5;
1890  }
1891  else
1892  {
1893  if (params->m_channelCondition->IsO2i ()) // outdoor to indoor
1894  {
1895  corrDis = 5;
1896  }
1897  else //LOS or NLOS
1898  {
1899  corrDis = 10;
1900  }
1901  }
1902  double R;
1903  if (m_blockerSpeed > 1e-6) // speed not equal to 0
1904  {
1905  double corrT = corrDis / m_blockerSpeed;
1906  R = exp (-1 * (deltaX / corrDis + (Now ().GetSeconds () - params->m_generatedTime.GetSeconds ()) / corrT));
1907  }
1908  else
1909  {
1910  R = exp (-1 * (deltaX / corrDis));
1911  }
1912 
1913  NS_LOG_INFO ("Distance change:" << deltaX << " Speed:" << m_blockerSpeed
1914  << " Time difference:" << Now ().GetSeconds () - params->m_generatedTime.GetSeconds ()
1915  << " correlation:" << R);
1916 
1917  //In order to generate correlated uniform random variables, we first generate correlated normal random variables and map the normal RV to uniform RV.
1918  //Notice the correlation will change if the RV is transformed from normal to uniform.
1919  //To compensate the distortion, the correlation of the normal RV is computed
1920  //such that the uniform RV would have the desired correlation when transformed from normal RV.
1921 
1922  //The following formula was obtained from MATLAB numerical simulation.
1923 
1924  if (R * R * (-0.069) + R * 1.074 - 0.002 < 1) //transform only when the correlation of normal RV is smaller than 1
1925  {
1926  R = R * R * (-0.069) + R * 1.074 - 0.002;
1927  }
1928  for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
1929  {
1930 
1931  //Generate a new correlated normal RV with the following formula
1932  params->m_nonSelfBlocking[blockInd][PHI_INDEX] =
1933  R * params->m_nonSelfBlocking[blockInd][PHI_INDEX] + sqrt (1 - R * R) * m_normalRv->GetValue ();
1934  }
1935  }
1936 
1937  }
1938 
1939  //step c: Determine the attenuation of each blocker due to blockers
1940  for (uint8_t cInd = 0; cInd < clusterNum; cInd++)
1941  {
1942  NS_ASSERT_MSG (clusterAOA[cInd] >= 0 && clusterAOA[cInd] <= 360, "the AOA should be the range of [0,360]");
1943  NS_ASSERT_MSG (clusterZOA[cInd] >= 0 && clusterZOA[cInd] <= 180, "the ZOA should be the range of [0,180]");
1944 
1945  //check self blocking
1946  NS_LOG_INFO ("AOA=" << clusterAOA[cInd] << " Block Region[" << phi_sb - x_sb / 2 << "," << phi_sb + x_sb / 2 << "]");
1947  NS_LOG_INFO ("ZOA=" << clusterZOA[cInd] << " Block Region[" << theta_sb - y_sb / 2 << "," << theta_sb + y_sb / 2 << "]");
1948  if ( std::abs (clusterAOA[cInd] - phi_sb) < (x_sb / 2) && std::abs (clusterZOA[cInd] - theta_sb) < (y_sb / 2))
1949  {
1950  powerAttenuation[cInd] += 30; //anttenuate by 30 dB.
1951  NS_LOG_INFO ("Cluster[" << (int)cInd << "] is blocked by self blocking region and reduce 30 dB power,"
1952  "the attenuation is [" << powerAttenuation[cInd] << " dB]");
1953  }
1954 
1955  //check non-self blocking
1956  double phiK, xK, thetaK, yK;
1957  for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
1958  {
1959  //The normal RV is transformed to uniform RV with the desired correlation.
1960  phiK = (0.5 * erfc (-1 * params->m_nonSelfBlocking[blockInd][PHI_INDEX] / sqrt (2))) * 360;
1961  while (phiK > 360)
1962  {
1963  phiK -= 360;
1964  }
1965 
1966  while (phiK < 0)
1967  {
1968  phiK += 360;
1969  }
1970 
1971  xK = params->m_nonSelfBlocking[blockInd][X_INDEX];
1972  thetaK = params->m_nonSelfBlocking[blockInd][THETA_INDEX];
1973  yK = params->m_nonSelfBlocking[blockInd][Y_INDEX];
1974  NS_LOG_INFO ("AOA=" << clusterAOA[cInd] << " Block Region[" << phiK - xK << "," << phiK + xK << "]");
1975  NS_LOG_INFO ("ZOA=" << clusterZOA[cInd] << " Block Region[" << thetaK - yK << "," << thetaK + yK << "]");
1976 
1977  if ( std::abs (clusterAOA[cInd] - phiK) < (xK)
1978  && std::abs (clusterZOA[cInd] - thetaK) < (yK))
1979  {
1980  double A1 = clusterAOA[cInd] - (phiK + xK / 2); //(7.6-24)
1981  double A2 = clusterAOA[cInd] - (phiK - xK / 2); //(7.6-25)
1982  double Z1 = clusterZOA[cInd] - (thetaK + yK / 2); //(7.6-26)
1983  double Z2 = clusterZOA[cInd] - (thetaK - yK / 2); //(7.6-27)
1984  int signA1, signA2, signZ1, signZ2;
1985  //draw sign for the above parameters according to table 7.6.4.1-3 Description of signs
1986  if (xK / 2 < clusterAOA[cInd] - phiK && clusterAOA[cInd] - phiK <= xK)
1987  {
1988  signA1 = -1;
1989  }
1990  else
1991  {
1992  signA1 = 1;
1993  }
1994  if (-1 * xK < clusterAOA[cInd] - phiK && clusterAOA[cInd] - phiK <= -1 * xK / 2)
1995  {
1996  signA2 = -1;
1997  }
1998  else
1999  {
2000  signA2 = 1;
2001  }
2002 
2003  if (yK / 2 < clusterZOA[cInd] - thetaK && clusterZOA[cInd] - thetaK <= yK)
2004  {
2005  signZ1 = -1;
2006  }
2007  else
2008  {
2009  signZ1 = 1;
2010  }
2011  if (-1 * yK < clusterZOA[cInd] - thetaK && clusterZOA[cInd] - thetaK <= -1 * yK / 2)
2012  {
2013  signZ2 = -1;
2014  }
2015  else
2016  {
2017  signZ2 = 1;
2018  }
2019  double lambda = 3e8 / m_frequency;
2020  double F_A1 = atan (signA1 * M_PI / 2 * sqrt (M_PI / lambda *
2021  params->m_nonSelfBlocking[blockInd][R_INDEX] * (1 / cos (A1 * M_PI / 180) - 1))) / M_PI; //(7.6-23)
2022  double F_A2 = atan (signA2 * M_PI / 2 * sqrt (M_PI / lambda *
2023  params->m_nonSelfBlocking[blockInd][R_INDEX] * (1 / cos (A2 * M_PI / 180) - 1))) / M_PI;
2024  double F_Z1 = atan (signZ1 * M_PI / 2 * sqrt (M_PI / lambda *
2025  params->m_nonSelfBlocking[blockInd][R_INDEX] * (1 / cos (Z1 * M_PI / 180) - 1))) / M_PI;
2026  double F_Z2 = atan (signZ2 * M_PI / 2 * sqrt (M_PI / lambda *
2027  params->m_nonSelfBlocking[blockInd][R_INDEX] * (1 / cos (Z2 * M_PI / 180) - 1))) / M_PI;
2028  double L_dB = -20 * log10 (1 - (F_A1 + F_A2) * (F_Z1 + F_Z2)); //(7.6-22)
2029  powerAttenuation[cInd] += L_dB;
2030  NS_LOG_INFO ("Cluster[" << (int)cInd << "] is blocked by no-self blocking, "
2031  "the loss is [" << L_dB << "]" << " dB");
2032 
2033  }
2034  }
2035  }
2036  return powerAttenuation;
2037 }
2038 
2039 
2040 void
2041 ThreeGppChannelModel::Shuffle (double * first, double * last) const
2042 {
2043  for (auto i = (last - first) - 1; i > 0; --i)
2044  {
2045  std::swap (first[i], first[m_uniformRvShuffle->GetInteger (0, i)]);
2046  }
2047 }
2048 
2049 int64_t
2050 ThreeGppChannelModel::AssignStreams (int64_t stream)
2051 {
2052  NS_LOG_FUNCTION (this << stream);
2053  m_normalRv->SetStream (stream);
2054  m_uniformRv->SetStream (stream + 1);
2055  m_uniformRvShuffle->SetStream (stream + 2);
2056  return 3;
2057 }
2058 
2059 } // namespace ns3
ns3::offSetAlpha
static const double offSetAlpha[20]
Definition: three-gpp-channel-model.cc:44
ns3::sqrtC_UMa_LOS
static const double sqrtC_UMa_LOS[7][7]
Definition: three-gpp-channel-model.cc:86
ns3::Time::As
TimeWithUnit As(const enum Unit unit=Time::AUTO) const
Attach a unit to a Time, to facilitate output in a specific unit.
Definition: time.cc:429
ns3::Time::NS
nanosecond
Definition: nstime.h:118
ns3::MatrixBasedChannelModel
This is an interface for a channel model that can be described by a channel matrix, e.g., the 3GPP Spatial Channel Models, which is generally used in combination with antenna arrays.
Definition: matrix-based-channel-model.h:43
ns3::ThreeGppChannelModel::m_channelMap
std::unordered_map< uint32_t, Ptr< ThreeGppChannelMatrix > > m_channelMap
map containing the channel realizations
Definition: three-gpp-channel-model.h:244
ns3::Ptr
Smart pointer class similar to boost::intrusive_ptr.
Definition: ptr.h:73
NS_LOG_FUNCTION
#define NS_LOG_FUNCTION(parameters)
If log level LOG_FUNCTION is enabled, this macro will output all input parameters separated by "...
Definition: log-macros-enabled.h:244
ns3::RandomVariableStream::SetStream
void SetStream(int64_t stream)
Specifies the stream number for the RngStream.
Definition: random-variable-stream.cc:100
ns3::BooleanValue
AttributeValue implementation for Boolean.
Definition: boolean.h:36
ns3::sqrtC_office_NLOS
static const double sqrtC_office_NLOS[6][6]
Definition: three-gpp-channel-model.cc:154
ns3::Node::GetId
uint32_t GetId(void) const
Definition: node.cc:109
ns3::UniformRandomVariable::GetInteger
uint32_t GetInteger(uint32_t min, uint32_t max)
Get the next random value, as an unsigned integer in the specified range .
Definition: random-variable-stream.cc:193
NS_OBJECT_ENSURE_REGISTERED
#define NS_OBJECT_ENSURE_REGISTERED(type)
Register an Object subclass with the TypeId system.
Definition: object-base.h:45
ns3::StringValue
Hold variables of type string.
Definition: string.h:41
min
#define min(a, b)
Definition: 80211b.c:42
ns3::sqrtC_UMi_O2I
static const double sqrtC_UMi_O2I[6][6]
Definition: three-gpp-channel-model.cc:135
ns3::ThreeGppChannelModel::DoDispose
void DoDispose() override
Destructor implementation.
Definition: three-gpp-channel-model.cc:180
ns3::ThreeGppChannelModel::THETA_INDEX
static const uint8_t THETA_INDEX
index of the THETA value in the m_nonSelfBlocking array
Definition: three-gpp-channel-model.h:261
ns3::MakeBooleanAccessor
Ptr< const AttributeAccessor > MakeBooleanAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: boolean.h:85
ns3::ThreeGppChannelModel::m_normalRv
Ptr< NormalRandomVariable > m_normalRv
normal random variable
Definition: three-gpp-channel-model.h:250
ns3::ThreeGppChannelModel::GetScenario
std::string GetScenario(void) const
Returns the propagation scenario.
Definition: three-gpp-channel-model.cc:285
ns3::ThreeGppChannelModel::m_uniformRvShuffle
Ptr< UniformRandomVariable > m_uniformRvShuffle
uniform random variable used to shuffle array in GetNewChannel
Definition: three-gpp-channel-model.h:251
NS_ASSERT
#define NS_ASSERT(condition)
At runtime, in debugging builds, if this condition is not true, the program prints the source file...
Definition: assert.h:67
ns3::IntegerValue
Hold a signed integer type.
Definition: integer.h:44
ns3::MatrixBasedChannelModel::Complex3DVector
std::vector< Complex2DVector > Complex3DVector
type definition for complex 3D matrices
Definition: matrix-based-channel-model.h:55
NS_LOG_COMPONENT_DEFINE
#define NS_LOG_COMPONENT_DEFINE(name)
Define a Log component with a specific name.
Definition: log.h:205
ns3::MilliSeconds
Time MilliSeconds(uint64_t value)
Construct a Time in the indicated unit.
Definition: nstime.h:1297
NS_LOG_INFO
#define NS_LOG_INFO(msg)
Use NS_LOG to output a message of level LOG_INFO.
Definition: log.h:281
NS_FATAL_ERROR
#define NS_FATAL_ERROR(msg)
Report a fatal error with a message and terminate.
Definition: fatal-error.h:165
ns3::Angles::theta
double theta
the inclination angle in radians
Definition: angles.h:117
ns3::MatrixBasedChannelModel::DoubleVector
std::vector< double > DoubleVector
type definition for vectors of doubles
Definition: matrix-based-channel-model.h:51
ns3::ThreeGppChannelModel::R_INDEX
static const uint8_t R_INDEX
index of the R value in the m_nonSelfBlocking array
Definition: three-gpp-channel-model.h:263
ns3::MatrixBasedChannelModel::ChannelMatrix::m_delay
DoubleVector m_delay
cluster delay in nanoseconds.
Definition: matrix-based-channel-model.h:64
ns3::ThreeGppChannelModel::SetChannelConditionModel
void SetChannelConditionModel(Ptr< ChannelConditionModel > model)
Set the channel condition model.
Definition: three-gpp-channel-model.cc:244
ns3::ThreeGppChannelModel::m_channelConditionModel
Ptr< ChannelConditionModel > m_channelConditionModel
the channel condition model
Definition: three-gpp-channel-model.h:248
ns3::ThreeGppChannelModel::m_scenario
std::string m_scenario
the 3GPP scenario
Definition: three-gpp-channel-model.h:247
ns3::sqrtC_UMa_O2I
static const double sqrtC_UMa_O2I[6][6]
Definition: three-gpp-channel-model.cc:106
ns3::MatrixBasedChannelModel::ChannelMatrix::m_nodeIds
std::pair< uint32_t, uint32_t > m_nodeIds
the first element is the s-node ID (the transmitter when the channel was generated), the second element is the u-node ID (the receiver when the channel was generated)
Definition: matrix-based-channel-model.h:67
ns3::ThreeGppChannelModel::X_INDEX
static const uint8_t X_INDEX
index of the X value in the m_nonSelfBlocking array
Definition: three-gpp-channel-model.h:260
ns3::MakeIntegerAccessor
Ptr< const AttributeAccessor > MakeIntegerAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: integer.h:45
ns3::MakeStringChecker
Ptr< const AttributeChecker > MakeStringChecker(void)
Definition: string.cc:30
ns3::sqrtC_office_LOS
static const double sqrtC_office_LOS[7][7]
Definition: three-gpp-channel-model.cc:144
ns3::ThreeGppChannelModel::SetFrequency
void SetFrequency(double f)
Sets the center frequency of the model.
Definition: three-gpp-channel-model.cc:258
ns3::ThreeGppChannelModel::GetFrequency
double GetFrequency(void) const
Returns the center frequency.
Definition: three-gpp-channel-model.cc:266
ns3::MakePointerAccessor
Ptr< const AttributeAccessor > MakePointerAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: pointer.h:227
max
#define max(a, b)
Definition: 80211b.c:43
ns3::Time::IsZero
bool IsZero(void) const
Exactly equivalent to t == 0.
Definition: nstime.h:301
ns3::TimeValue
AttributeValue implementation for Time.
Definition: nstime.h:1353
ns3::ThreeGppChannelModel::m_numNonSelfBlocking
uint16_t m_numNonSelfBlocking
number of non-self-blocking regions
Definition: three-gpp-channel-model.h:255
ns3::sqrtC_RMa_LOS
static const double sqrtC_RMa_LOS[7][7]
Definition: three-gpp-channel-model.cc:58
ns3::ThreeGppChannelModel::m_updatePeriod
Time m_updatePeriod
the channel update period
Definition: three-gpp-channel-model.h:245
ns3::ThreeGppChannelModel::AssignStreams
int64_t AssignStreams(int64_t stream)
Assign a fixed random variable stream number to the random variables used by this model...
Definition: three-gpp-channel-model.cc:2050
ns3::ThreeGppChannelModel::GetChannel
Ptr< const ChannelMatrix > GetChannel(Ptr< const MobilityModel > aMob, Ptr< const MobilityModel > bMob, Ptr< const ThreeGppAntennaArrayModel > aAntenna, Ptr< const ThreeGppAntennaArrayModel > bAntenna) override
Looks for the channel matrix associated to the aMob and bMob pair in m_channelMap.
Definition: three-gpp-channel-model.cc:971
ns3::ThreeGppChannelModel::m_portraitMode
bool m_portraitMode
true if potrait mode, false if landscape
Definition: three-gpp-channel-model.h:256
ns3::ThreeGppChannelModel::ThreeGppChannelMatrix::m_channelCondition
Ptr< const ChannelCondition > m_channelCondition
the channel condition
Definition: three-gpp-channel-model.h:144
ns3::MatrixBasedChannelModel::ChannelMatrix::m_channel
Complex3DVector m_channel
channel matrix H[u][s][n].
Definition: matrix-based-channel-model.h:63
ns3::ThreeGppChannelModel::ThreeGppChannelMatrix::m_clusterPhase
MatrixBasedChannelModel::Double3DVector m_clusterPhase
the initial random phases
Definition: three-gpp-channel-model.h:156
ns3::ThreeGppChannelModel::m_blockerSpeed
double m_blockerSpeed
the blocker speed
Definition: three-gpp-channel-model.h:257
f
double f(double x, void *params)
Definition: 80211b.c:70
ns3::ThreeGppChannelModel::GetThreeGppTable
virtual Ptr< const ParamsTable > GetThreeGppTable(Ptr< const ChannelCondition > channelCondition, double hBS, double hUT, double distance2D) const
Get the parameters needed to apply the channel generation procedure.
Definition: three-gpp-channel-model.cc:292
ns3
Every class exported by the ns3 library is enclosed in the ns3 namespace.
ns3::PointerValue
Hold objects of type Ptr<T>.
Definition: pointer.h:36
ns3::UniformRandomVariable::GetValue
double GetValue(double min, double max)
Get the next random value, as a double in the specified range .
Definition: random-variable-stream.cc:182
ns3::ThreeGppChannelModel::m_uniformRv
Ptr< UniformRandomVariable > m_uniformRv
uniform random variable
Definition: three-gpp-channel-model.h:249
ns3::MatrixBasedChannelModel::GetKey
static constexpr uint32_t GetKey(uint32_t x1, uint32_t x2)
Calculate the channel key using the Cantor function.
Definition: matrix-based-channel-model.h:120
ns3::MatrixBasedChannelModel::Double2DVector
std::vector< DoubleVector > Double2DVector
type definition for matrices of doubles
Definition: matrix-based-channel-model.h:52
ns3::MakeTimeAccessor
Ptr< const AttributeAccessor > MakeTimeAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: nstime.h:1354
ns3::Simulator::Now
static Time Now(void)
Return the current simulation virtual time.
Definition: simulator.cc:195
ns3::sqrtC_RMa_NLOS
static const double sqrtC_RMa_NLOS[6][6]
Definition: three-gpp-channel-model.cc:68
ns3::MatrixBasedChannelModel::ChannelMatrix::m_angle
Double2DVector m_angle
cluster angle angle[direction][n], where direction = 0(AOA), 1(ZOA), 2(AOD), 3(ZOD) in degree...
Definition: matrix-based-channel-model.h:65
ns3::sqrtC_RMa_O2I
static const double sqrtC_RMa_O2I[6][6]
Definition: three-gpp-channel-model.cc:77
ns3::MatrixBasedChannelModel::Double3DVector
std::vector< Double2DVector > Double3DVector
type definition for 3D matrices of doubles
Definition: matrix-based-channel-model.h:53
NS_ASSERT_MSG
#define NS_ASSERT_MSG(condition, message)
At runtime, in debugging builds, if this condition is not true, the program prints the message to out...
Definition: assert.h:88
ns3::MakeDoubleAccessor
Ptr< const AttributeAccessor > MakeDoubleAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: double.h:42
ns3::max
double max(double x, double y)
Definition: cobalt-queue-disc.cc:122
ns3::ThreeGppChannelModel::SetScenario
void SetScenario(const std::string &scenario)
Sets the propagation scenario.
Definition: three-gpp-channel-model.cc:273
ns3::ThreeGppChannelModel::PHI_INDEX
static const uint8_t PHI_INDEX
index of the PHI value in the m_nonSelfBlocking array
Definition: three-gpp-channel-model.h:259
ns3::ThreeGppChannelModel::m_blockage
bool m_blockage
enables the blockage model A
Definition: three-gpp-channel-model.h:254
NS_LOG_WARN
#define NS_LOG_WARN(msg)
Use NS_LOG to output a message of level LOG_WARN.
Definition: log.h:265
ns3::ThreeGppChannelModel::GetNewChannel
Ptr< ThreeGppChannelMatrix > GetNewChannel(Vector locUT, Ptr< const ChannelCondition > channelCondition, Ptr< const ThreeGppAntennaArrayModel > sAntenna, Ptr< const ThreeGppAntennaArrayModel > uAntenna, Angles &uAngle, Angles &sAngle, double dis2D, double hBS, double hUT) const
Compute the channel matrix between two devices using the procedure described in 3GPP TR 38...
Definition: three-gpp-channel-model.cc:1040
ns3::ThreeGppChannelModel::GetTypeId
static TypeId GetTypeId()
Get the type ID.
Definition: three-gpp-channel-model.cc:188
ns3::ThreeGppChannelModel::ChannelMatrixNeedsUpdate
bool ChannelMatrixNeedsUpdate(Ptr< const ThreeGppChannelMatrix > channelMatrix, Ptr< const ChannelCondition > channelCondition) const
Check if the channel matrix has to be updated.
Definition: three-gpp-channel-model.cc:947
ns3::Node
A network Node.
Definition: node.h:56
NS_LOG_DEBUG
#define NS_LOG_DEBUG(msg)
Use NS_LOG to output a message of level LOG_DEBUG.
Definition: log.h:273
ns3::MakeBooleanChecker
Ptr< const AttributeChecker > MakeBooleanChecker(void)
Definition: boolean.cc:121
ns3::Angles::phi
double phi
the azimuth angle in radians
Definition: angles.h:111
ns3::ThreeGppChannelModel::Shuffle
void Shuffle(double *first, double *last) const
Shuffle the elements of a simple sequence container of type double.
Definition: three-gpp-channel-model.cc:2041
ns3::min
double min(double x, double y)
Definition: cobalt-queue-disc.cc:117
ns3::MakeTimeChecker
Ptr< const AttributeChecker > MakeTimeChecker(const Time min, const Time max)
Helper to make a Time checker with bounded range.
Definition: time.cc:533
ns3::Now
Time Now(void)
create an ns3::Time instance which contains the current simulation time.
Definition: simulator.cc:287
ns3::Angles
struct holding the azimuth and inclination angles of spherical coordinates.
Definition: angles.h:71
ns3::Object
A base class which provides memory management and object aggregation.
Definition: object.h:87
ns3::MakeStringAccessor
Ptr< const AttributeAccessor > MakeStringAccessor(T1 a1)
Create an AttributeAccessor for a class data member, or a lone class get functor or set method...
Definition: string.h:42
ns3::sqrtC_UMa_NLOS
static const double sqrtC_UMa_NLOS[6][6]
Definition: three-gpp-channel-model.cc:97
first
Definition: first.py:1
ns3::DoubleValue
This class can be used to hold variables of floating point type such as &#39;double&#39; or &#39;float&#39;...
Definition: double.h:41
ns3::ThreeGppChannelModel::GetChannelConditionModel
Ptr< ChannelConditionModel > GetChannelConditionModel() const
Get the associated channel condition model.
Definition: three-gpp-channel-model.cc:251
ns3::ThreeGppChannelModel::Y_INDEX
static const uint8_t Y_INDEX
index of the Y value in the m_nonSelfBlocking array
Definition: three-gpp-channel-model.h:262
ns3::TypeId
a unique identifier for an interface.
Definition: type-id.h:58
ns3::sqrtC_UMi_LOS
static const double sqrtC_UMi_LOS[7][7]
Definition: three-gpp-channel-model.cc:116
ns3::ThreeGppChannelModel::CalcAttenuationOfBlockage
DoubleVector CalcAttenuationOfBlockage(Ptr< ThreeGppChannelMatrix > params, const DoubleVector &clusterAOA, const DoubleVector &clusterZOA) const
Applies the blockage model A described in 3GPP TR 38.901.
Definition: three-gpp-channel-model.cc:1820
ns3::TypeId::SetParent
TypeId SetParent(TypeId tid)
Set the parent TypeId.
Definition: type-id.cc:923
ns3::sqrtC_UMi_NLOS
static const double sqrtC_UMi_NLOS[6][6]
Definition: three-gpp-channel-model.cc:126
ns3::ThreeGppChannelModel::m_frequency
double m_frequency
the operating frequency
Definition: three-gpp-channel-model.h:246