A Discrete-Event Network Simulator
API
Loading...
Searching...
No Matches
three-gpp-channel-model.cc
Go to the documentation of this file.
1/*
2 * Copyright (c) 2019 SIGNET Lab, Department of Information Engineering,
3 * University of Padova
4 * Copyright (c) 2015, NYU WIRELESS, Tandon School of Engineering,
5 * New York University
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License version 2 as
9 * published by the Free Software Foundation;
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 *
20 */
21
23
24#include "ns3/double.h"
25#include "ns3/integer.h"
26#include "ns3/log.h"
27#include "ns3/mobility-model.h"
28#include "ns3/node.h"
29#include "ns3/phased-array-model.h"
30#include "ns3/pointer.h"
31#include "ns3/string.h"
32#include <ns3/simulator.h>
33
34#include <algorithm>
35#include <random>
36
37namespace ns3
38{
39
40NS_LOG_COMPONENT_DEFINE("ThreeGppChannelModel");
41
42NS_OBJECT_ENSURE_REGISTERED(ThreeGppChannelModel);
43
46static const double offSetAlpha[20] = {
47 0.0447, -0.0447, 0.1413, -0.1413, 0.2492, -0.2492, 0.3715, -0.3715, 0.5129, -0.5129,
48 0.6797, -0.6797, 0.8844, -0.8844, 1.1481, -1.1481, 1.5195, -1.5195, 2.1551, -2.1551,
49};
50
59static const double sqrtC_RMa_LOS[7][7] = {
60 {1, 0, 0, 0, 0, 0, 0},
61 {0, 1, 0, 0, 0, 0, 0},
62 {-0.5, 0, 0.866025, 0, 0, 0, 0},
63 {0, 0, 0, 1, 0, 0, 0},
64 {0, 0, 0, 0, 1, 0, 0},
65 {0.01, 0, -0.0519615, 0.73, -0.2, 0.651383, 0},
66 {-0.17, -0.02, 0.21362, -0.14, 0.24, 0.142773, 0.909661},
67};
68
78static const double sqrtC_RMa_NLOS[6][6] = {
79 {1, 0, 0, 0, 0, 0},
80 {-0.5, 0.866025, 0, 0, 0, 0},
81 {0.6, -0.11547, 0.791623, 0, 0, 0},
82 {0, 0, 0, 1, 0, 0},
83 {-0.04, -0.138564, 0.540662, -0.18, 0.809003, 0},
84 {-0.25, -0.606218, -0.240013, 0.26, -0.231685, 0.625392},
85};
86
95static const double sqrtC_RMa_O2I[6][6] = {
96 {1, 0, 0, 0, 0, 0},
97 {0, 1, 0, 0, 0, 0},
98 {0, 0, 1, 0, 0, 0},
99 {0, 0, -0.7, 0.714143, 0, 0},
100 {0, 0, 0.66, -0.123225, 0.741091, 0},
101 {0, 0, 0.47, 0.152631, -0.393194, 0.775373},
102};
103
112static const double sqrtC_UMa_LOS[7][7] = {
113 {1, 0, 0, 0, 0, 0, 0},
114 {0, 1, 0, 0, 0, 0, 0},
115 {-0.4, -0.4, 0.824621, 0, 0, 0, 0},
116 {-0.5, 0, 0.242536, 0.83137, 0, 0, 0},
117 {-0.5, -0.2, 0.630593, -0.484671, 0.278293, 0, 0},
118 {0, 0, -0.242536, 0.672172, 0.642214, 0.27735, 0},
119 {-0.8, 0, -0.388057, -0.367926, 0.238537, -3.58949e-15, 0.130931},
120};
121
131static const double sqrtC_UMa_NLOS[6][6] = {
132 {1, 0, 0, 0, 0, 0},
133 {-0.4, 0.916515, 0, 0, 0, 0},
134 {-0.6, 0.174574, 0.78072, 0, 0, 0},
135 {0, 0.654654, 0.365963, 0.661438, 0, 0},
136 {0, -0.545545, 0.762422, 0.118114, 0.327327, 0},
137 {-0.4, -0.174574, -0.396459, 0.392138, 0.49099, 0.507445},
138};
139
148static const double sqrtC_UMa_O2I[6][6] = {
149 {1, 0, 0, 0, 0, 0},
150 {-0.5, 0.866025, 0, 0, 0, 0},
151 {0.2, 0.57735, 0.791623, 0, 0, 0},
152 {0, 0.46188, -0.336861, 0.820482, 0, 0},
153 {0, -0.69282, 0.252646, 0.493742, 0.460857, 0},
154 {0, -0.23094, 0.16843, 0.808554, -0.220827, 0.464515},
155
156};
157
166static const double sqrtC_UMi_LOS[7][7] = {
167 {1, 0, 0, 0, 0, 0, 0},
168 {0.5, 0.866025, 0, 0, 0, 0, 0},
169 {-0.4, -0.57735, 0.711805, 0, 0, 0, 0},
170 {-0.5, 0.057735, 0.468293, 0.726201, 0, 0, 0},
171 {-0.4, -0.11547, 0.805464, -0.23482, 0.350363, 0, 0},
172 {0, 0, 0, 0.688514, 0.461454, 0.559471, 0},
173 {0, 0, 0.280976, 0.231921, -0.490509, 0.11916, 0.782603},
174};
175
185static const double sqrtC_UMi_NLOS[6][6] = {
186 {1, 0, 0, 0, 0, 0},
187 {-0.7, 0.714143, 0, 0, 0, 0},
188 {0, 0, 1, 0, 0, 0},
189 {-0.4, 0.168034, 0, 0.90098, 0, 0},
190 {0, -0.70014, 0.5, 0.130577, 0.4927, 0},
191 {0, 0, 0.5, 0.221981, -0.566238, 0.616522},
192};
193
202static const double sqrtC_UMi_O2I[6][6] = {
203 {1, 0, 0, 0, 0, 0},
204 {-0.5, 0.866025, 0, 0, 0, 0},
205 {0.2, 0.57735, 0.791623, 0, 0, 0},
206 {0, 0.46188, -0.336861, 0.820482, 0, 0},
207 {0, -0.69282, 0.252646, 0.493742, 0.460857, 0},
208 {0, -0.23094, 0.16843, 0.808554, -0.220827, 0.464515},
209};
210
219static const double sqrtC_office_LOS[7][7] = {
220 {1, 0, 0, 0, 0, 0, 0},
221 {0.5, 0.866025, 0, 0, 0, 0, 0},
222 {-0.8, -0.11547, 0.588784, 0, 0, 0, 0},
223 {-0.4, 0.23094, 0.520847, 0.717903, 0, 0, 0},
224 {-0.5, 0.288675, 0.73598, -0.348236, 0.0610847, 0, 0},
225 {0.2, -0.11547, 0.418943, 0.541106, 0.219905, 0.655744, 0},
226 {0.3, -0.057735, 0.73598, -0.348236, 0.0610847, -0.304997, 0.383375},
227};
228
238static const double sqrtC_office_NLOS[6][6] = {
239 {1, 0, 0, 0, 0, 0},
240 {-0.5, 0.866025, 0, 0, 0, 0},
241 {0, 0.46188, 0.886942, 0, 0, 0},
242 {-0.4, -0.23094, 0.120263, 0.878751, 0, 0},
243 {0, -0.311769, 0.55697, -0.249198, 0.728344, 0},
244 {0, -0.069282, 0.295397, 0.430696, 0.468462, 0.709214},
245};
246
248{
249 NS_LOG_FUNCTION(this);
250 m_uniformRv = CreateObject<UniformRandomVariable>();
251 m_uniformRvShuffle = CreateObject<UniformRandomVariable>();
252 m_uniformRvDoppler = CreateObject<UniformRandomVariable>();
253
254 m_normalRv = CreateObject<NormalRandomVariable>();
255 m_normalRv->SetAttribute("Mean", DoubleValue(0.0));
256 m_normalRv->SetAttribute("Variance", DoubleValue(1.0));
257}
258
260{
261 NS_LOG_FUNCTION(this);
262}
263
264void
266{
267 NS_LOG_FUNCTION(this);
269 {
270 m_channelConditionModel->Dispose();
271 }
272 m_channelMatrixMap.clear();
273 m_channelParamsMap.clear();
274 m_channelConditionModel = nullptr;
275}
276
277TypeId
279{
280 static TypeId tid =
281 TypeId("ns3::ThreeGppChannelModel")
282 .SetGroupName("Spectrum")
284 .AddConstructor<ThreeGppChannelModel>()
285 .AddAttribute("Frequency",
286 "The operating Frequency in Hz",
287 DoubleValue(500.0e6),
290 MakeDoubleChecker<double>())
291 .AddAttribute(
292 "Scenario",
293 "The 3GPP scenario (RMa, UMa, UMi-StreetCanyon, InH-OfficeOpen, InH-OfficeMixed)",
294 StringValue("UMa"),
298 .AddAttribute("ChannelConditionModel",
299 "Pointer to the channel condition model",
300 PointerValue(),
303 MakePointerChecker<ChannelConditionModel>())
304 .AddAttribute("UpdatePeriod",
305 "Specify the channel coherence time",
309 // attributes for the blockage model
310 .AddAttribute("Blockage",
311 "Enable blockage model A (sec 7.6.4.1)",
312 BooleanValue(false),
315 .AddAttribute("NumNonselfBlocking",
316 "number of non-self-blocking regions",
317 IntegerValue(4),
319 MakeIntegerChecker<uint16_t>())
320 .AddAttribute("PortraitMode",
321 "true for portrait mode, false for landscape mode",
322 BooleanValue(true),
325 .AddAttribute("BlockerSpeed",
326 "The speed of moving blockers, the unit is m/s",
327 DoubleValue(1),
329 MakeDoubleChecker<double>())
330 .AddAttribute("vScatt",
331 "Maximum speed of the vehicle in the layout (see 3GPP TR 37.885 v15.3.0, "
332 "Sec. 6.2.3)."
333 "Used to compute the additional contribution for the Doppler of"
334 "delayed (reflected) paths",
335 DoubleValue(0.0),
337 MakeDoubleChecker<double>(0.0))
338
339 ;
340 return tid;
341}
342
343void
345{
346 NS_LOG_FUNCTION(this);
348}
349
352{
353 NS_LOG_FUNCTION(this);
355}
356
357void
359{
360 NS_LOG_FUNCTION(this);
361 NS_ASSERT_MSG(f >= 500.0e6 && f <= 100.0e9,
362 "Frequency should be between 0.5 and 100 GHz but is " << f);
363 m_frequency = f;
364}
365
366double
368{
369 NS_LOG_FUNCTION(this);
370 return m_frequency;
371}
372
373void
374ThreeGppChannelModel::SetScenario(const std::string& scenario)
375{
376 NS_LOG_FUNCTION(this);
377 NS_ASSERT_MSG(scenario == "RMa" || scenario == "UMa" || scenario == "UMi-StreetCanyon" ||
378 scenario == "InH-OfficeOpen" || scenario == "InH-OfficeMixed" ||
379 scenario == "V2V-Urban" || scenario == "V2V-Highway",
380 "Unknown scenario, choose between: RMa, UMa, UMi-StreetCanyon, "
381 "InH-OfficeOpen, InH-OfficeMixed, V2V-Urban or V2V-Highway");
382 m_scenario = scenario;
383}
384
385std::string
387{
388 NS_LOG_FUNCTION(this);
389 return m_scenario;
390}
391
394 double hBS,
395 double hUT,
396 double distance2D) const
397{
398 NS_LOG_FUNCTION(this);
399
400 double fcGHz = m_frequency / 1.0e9;
401 Ptr<ParamsTable> table3gpp = Create<ParamsTable>();
402 // table3gpp includes the following parameters:
403 // numOfCluster, raysPerCluster, uLgDS, sigLgDS, uLgASD, sigLgASD,
404 // uLgASA, sigLgASA, uLgZSA, sigLgZSA, uLgZSD, sigLgZSD, offsetZOD,
405 // cDS, cASD, cASA, cZSA, uK, sigK, rTau, uXpr, sigXpr, shadowingStd
406
407 bool los = channelCondition->IsLos();
408 bool o2i = channelCondition->IsO2i();
409
410 // In NLOS case, parameter uK and sigK are not used and they are set to 0
411 if (m_scenario == "RMa")
412 {
413 if (los && !o2i)
414 {
415 // 3GPP mentioned that 3.91 ns should be used when the Cluster DS (cDS)
416 // entry is N/A.
417 table3gpp->m_numOfCluster = 11;
418 table3gpp->m_raysPerCluster = 20;
419 table3gpp->m_uLgDS = -7.49;
420 table3gpp->m_sigLgDS = 0.55;
421 table3gpp->m_uLgASD = 0.90;
422 table3gpp->m_sigLgASD = 0.38;
423 table3gpp->m_uLgASA = 1.52;
424 table3gpp->m_sigLgASA = 0.24;
425 table3gpp->m_uLgZSA = 0.47;
426 table3gpp->m_sigLgZSA = 0.40;
427 table3gpp->m_uLgZSD = 0.34;
428 table3gpp->m_sigLgZSD =
429 std::max(-1.0, -0.17 * (distance2D / 1000.0) - 0.01 * (hUT - 1.5) + 0.22);
430 table3gpp->m_offsetZOD = 0;
431 table3gpp->m_cDS = 3.91e-9;
432 table3gpp->m_cASD = 2;
433 table3gpp->m_cASA = 3;
434 table3gpp->m_cZSA = 3;
435 table3gpp->m_uK = 7;
436 table3gpp->m_sigK = 4;
437 table3gpp->m_rTau = 3.8;
438 table3gpp->m_uXpr = 12;
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_RMa_LOS[row][column];
447 }
448 }
449 }
450 else if (!los && !o2i)
451 {
452 table3gpp->m_numOfCluster = 10;
453 table3gpp->m_raysPerCluster = 20;
454 table3gpp->m_uLgDS = -7.43;
455 table3gpp->m_sigLgDS = 0.48;
456 table3gpp->m_uLgASD = 0.95;
457 table3gpp->m_sigLgASD = 0.45;
458 table3gpp->m_uLgASA = 1.52;
459 table3gpp->m_sigLgASA = 0.13;
460 table3gpp->m_uLgZSA = 0.58;
461 table3gpp->m_sigLgZSA = 0.37;
462 table3gpp->m_uLgZSD =
463 std::max(-1.0, -0.19 * (distance2D / 1000.0) - 0.01 * (hUT - 1.5) + 0.28);
464 table3gpp->m_sigLgZSD = 0.30;
465 table3gpp->m_offsetZOD = atan((35 - 3.5) / distance2D) - atan((35 - 1.5) / distance2D);
466 table3gpp->m_cDS = 3.91e-9;
467 table3gpp->m_cASD = 2;
468 table3gpp->m_cASA = 3;
469 table3gpp->m_cZSA = 3;
470 table3gpp->m_uK = 0;
471 table3gpp->m_sigK = 0;
472 table3gpp->m_rTau = 1.7;
473 table3gpp->m_uXpr = 7;
474 table3gpp->m_sigXpr = 3;
475 table3gpp->m_perClusterShadowingStd = 3;
476
477 for (uint8_t row = 0; row < 6; row++)
478 {
479 for (uint8_t column = 0; column < 6; column++)
480 {
481 table3gpp->m_sqrtC[row][column] = sqrtC_RMa_NLOS[row][column];
482 }
483 }
484 }
485 else // o2i
486 {
487 table3gpp->m_numOfCluster = 10;
488 table3gpp->m_raysPerCluster = 20;
489 table3gpp->m_uLgDS = -7.47;
490 table3gpp->m_sigLgDS = 0.24;
491 table3gpp->m_uLgASD = 0.67;
492 table3gpp->m_sigLgASD = 0.18;
493 table3gpp->m_uLgASA = 1.66;
494 table3gpp->m_sigLgASA = 0.21;
495 table3gpp->m_uLgZSA = 0.93;
496 table3gpp->m_sigLgZSA = 0.22;
497 table3gpp->m_uLgZSD =
498 std::max(-1.0, -0.19 * (distance2D / 1000.0) - 0.01 * (hUT - 1.5) + 0.28);
499 table3gpp->m_sigLgZSD = 0.30;
500 table3gpp->m_offsetZOD = atan((35 - 3.5) / distance2D) - atan((35 - 1.5) / distance2D);
501 table3gpp->m_cDS = 3.91e-9;
502 table3gpp->m_cASD = 2;
503 table3gpp->m_cASA = 3;
504 table3gpp->m_cZSA = 3;
505 table3gpp->m_uK = 0;
506 table3gpp->m_sigK = 0;
507 table3gpp->m_rTau = 1.7;
508 table3gpp->m_uXpr = 7;
509 table3gpp->m_sigXpr = 3;
510 table3gpp->m_perClusterShadowingStd = 3;
511
512 for (uint8_t row = 0; row < 6; row++)
513 {
514 for (uint8_t column = 0; column < 6; column++)
515 {
516 table3gpp->m_sqrtC[row][column] = sqrtC_RMa_O2I[row][column];
517 }
518 }
519 }
520 }
521 else if (m_scenario == "UMa")
522 {
523 if (los && !o2i)
524 {
525 table3gpp->m_numOfCluster = 12;
526 table3gpp->m_raysPerCluster = 20;
527 table3gpp->m_uLgDS = -6.955 - 0.0963 * log10(fcGHz);
528 table3gpp->m_sigLgDS = 0.66;
529 table3gpp->m_uLgASD = 1.06 + 0.1114 * log10(fcGHz);
530 table3gpp->m_sigLgASD = 0.28;
531 table3gpp->m_uLgASA = 1.81;
532 table3gpp->m_sigLgASA = 0.20;
533 table3gpp->m_uLgZSA = 0.95;
534 table3gpp->m_sigLgZSA = 0.16;
535 table3gpp->m_uLgZSD =
536 std::max(-0.5, -2.1 * distance2D / 1000.0 - 0.01 * (hUT - 1.5) + 0.75);
537 table3gpp->m_sigLgZSD = 0.40;
538 table3gpp->m_offsetZOD = 0;
539 table3gpp->m_cDS = std::max(0.25, -3.4084 * log10(fcGHz) + 6.5622) * 1e-9;
540 table3gpp->m_cASD = 5;
541 table3gpp->m_cASA = 11;
542 table3gpp->m_cZSA = 7;
543 table3gpp->m_uK = 9;
544 table3gpp->m_sigK = 3.5;
545 table3gpp->m_rTau = 2.5;
546 table3gpp->m_uXpr = 8;
547 table3gpp->m_sigXpr = 4;
548 table3gpp->m_perClusterShadowingStd = 3;
549
550 for (uint8_t row = 0; row < 7; row++)
551 {
552 for (uint8_t column = 0; column < 7; column++)
553 {
554 table3gpp->m_sqrtC[row][column] = sqrtC_UMa_LOS[row][column];
555 }
556 }
557 }
558 else
559 {
560 double uLgZSD = std::max(-0.5, -2.1 * distance2D / 1000.0 - 0.01 * (hUT - 1.5) + 0.9);
561
562 double afc = 0.208 * log10(fcGHz) - 0.782;
563 double bfc = 25;
564 double cfc = -0.13 * log10(fcGHz) + 2.03;
565 double efc = 7.66 * log10(fcGHz) - 5.96;
566
567 double offsetZOD = efc - std::pow(10, afc * log10(std::max(bfc, distance2D)) + cfc);
568
569 if (!los && !o2i)
570 {
571 table3gpp->m_numOfCluster = 20;
572 table3gpp->m_raysPerCluster = 20;
573 table3gpp->m_uLgDS = -6.28 - 0.204 * log10(fcGHz);
574 table3gpp->m_sigLgDS = 0.39;
575 table3gpp->m_uLgASD = 1.5 - 0.1144 * log10(fcGHz);
576 table3gpp->m_sigLgASD = 0.28;
577 table3gpp->m_uLgASA = 2.08 - 0.27 * log10(fcGHz);
578 table3gpp->m_sigLgASA = 0.11;
579 table3gpp->m_uLgZSA = -0.3236 * log10(fcGHz) + 1.512;
580 table3gpp->m_sigLgZSA = 0.16;
581 table3gpp->m_uLgZSD = uLgZSD;
582 table3gpp->m_sigLgZSD = 0.49;
583 table3gpp->m_offsetZOD = offsetZOD;
584 table3gpp->m_cDS = std::max(0.25, -3.4084 * log10(fcGHz) + 6.5622) * 1e-9;
585 table3gpp->m_cASD = 2;
586 table3gpp->m_cASA = 15;
587 table3gpp->m_cZSA = 7;
588 table3gpp->m_uK = 0;
589 table3gpp->m_sigK = 0;
590 table3gpp->m_rTau = 2.3;
591 table3gpp->m_uXpr = 7;
592 table3gpp->m_sigXpr = 3;
593 table3gpp->m_perClusterShadowingStd = 3;
594
595 for (uint8_t row = 0; row < 6; row++)
596 {
597 for (uint8_t column = 0; column < 6; column++)
598 {
599 table3gpp->m_sqrtC[row][column] = sqrtC_UMa_NLOS[row][column];
600 }
601 }
602 }
603 else //(o2i)
604 {
605 table3gpp->m_numOfCluster = 12;
606 table3gpp->m_raysPerCluster = 20;
607 table3gpp->m_uLgDS = -6.62;
608 table3gpp->m_sigLgDS = 0.32;
609 table3gpp->m_uLgASD = 1.25;
610 table3gpp->m_sigLgASD = 0.42;
611 table3gpp->m_uLgASA = 1.76;
612 table3gpp->m_sigLgASA = 0.16;
613 table3gpp->m_uLgZSA = 1.01;
614 table3gpp->m_sigLgZSA = 0.43;
615 table3gpp->m_uLgZSD = uLgZSD;
616 table3gpp->m_sigLgZSD = 0.49;
617 table3gpp->m_offsetZOD = offsetZOD;
618 table3gpp->m_cDS = 11e-9;
619 table3gpp->m_cASD = 5;
620 table3gpp->m_cASA = 8;
621 table3gpp->m_cZSA = 3;
622 table3gpp->m_uK = 0;
623 table3gpp->m_sigK = 0;
624 table3gpp->m_rTau = 2.2;
625 table3gpp->m_uXpr = 9;
626 table3gpp->m_sigXpr = 5;
627 table3gpp->m_perClusterShadowingStd = 4;
628
629 for (uint8_t row = 0; row < 6; row++)
630 {
631 for (uint8_t column = 0; column < 6; column++)
632 {
633 table3gpp->m_sqrtC[row][column] = sqrtC_UMa_O2I[row][column];
634 }
635 }
636 }
637 }
638 }
639 else if (m_scenario == "UMi-StreetCanyon")
640 {
641 if (los && !o2i)
642 {
643 table3gpp->m_numOfCluster = 12;
644 table3gpp->m_raysPerCluster = 20;
645 table3gpp->m_uLgDS = -0.24 * log10(1 + fcGHz) - 7.14;
646 table3gpp->m_sigLgDS = 0.38;
647 table3gpp->m_uLgASD = -0.05 * log10(1 + fcGHz) + 1.21;
648 table3gpp->m_sigLgASD = 0.41;
649 table3gpp->m_uLgASA = -0.08 * log10(1 + fcGHz) + 1.73;
650 table3gpp->m_sigLgASA = 0.014 * log10(1 + fcGHz) + 0.28;
651 table3gpp->m_uLgZSA = -0.1 * log10(1 + fcGHz) + 0.73;
652 table3gpp->m_sigLgZSA = -0.04 * log10(1 + fcGHz) + 0.34;
653 table3gpp->m_uLgZSD =
654 std::max(-0.21, -14.8 * distance2D / 1000.0 + 0.01 * std::abs(hUT - hBS) + 0.83);
655 table3gpp->m_sigLgZSD = 0.35;
656 table3gpp->m_offsetZOD = 0;
657 table3gpp->m_cDS = 5e-9;
658 table3gpp->m_cASD = 3;
659 table3gpp->m_cASA = 17;
660 table3gpp->m_cZSA = 7;
661 table3gpp->m_uK = 9;
662 table3gpp->m_sigK = 5;
663 table3gpp->m_rTau = 3;
664 table3gpp->m_uXpr = 9;
665 table3gpp->m_sigXpr = 3;
666 table3gpp->m_perClusterShadowingStd = 3;
667
668 for (uint8_t row = 0; row < 7; row++)
669 {
670 for (uint8_t column = 0; column < 7; column++)
671 {
672 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
673 }
674 }
675 }
676 else
677 {
678 double uLgZSD =
679 std::max(-0.5, -3.1 * distance2D / 1000.0 + 0.01 * std::max(hUT - hBS, 0.0) + 0.2);
680 double offsetZOD = -1 * std::pow(10, -1.5 * log10(std::max(10.0, distance2D)) + 3.3);
681 if (!los && !o2i)
682 {
683 table3gpp->m_numOfCluster = 19;
684 table3gpp->m_raysPerCluster = 20;
685 table3gpp->m_uLgDS = -0.24 * log10(1 + fcGHz) - 6.83;
686 table3gpp->m_sigLgDS = 0.16 * log10(1 + fcGHz) + 0.28;
687 table3gpp->m_uLgASD = -0.23 * log10(1 + fcGHz) + 1.53;
688 table3gpp->m_sigLgASD = 0.11 * log10(1 + fcGHz) + 0.33;
689 table3gpp->m_uLgASA = -0.08 * log10(1 + fcGHz) + 1.81;
690 table3gpp->m_sigLgASA = 0.05 * log10(1 + fcGHz) + 0.3;
691 table3gpp->m_uLgZSA = -0.04 * log10(1 + fcGHz) + 0.92;
692 table3gpp->m_sigLgZSA = -0.07 * log10(1 + fcGHz) + 0.41;
693 table3gpp->m_uLgZSD = uLgZSD;
694 table3gpp->m_sigLgZSD = 0.35;
695 table3gpp->m_offsetZOD = offsetZOD;
696 table3gpp->m_cDS = 11e-9;
697 table3gpp->m_cASD = 10;
698 table3gpp->m_cASA = 22;
699 table3gpp->m_cZSA = 7;
700 table3gpp->m_uK = 0;
701 table3gpp->m_sigK = 0;
702 table3gpp->m_rTau = 2.1;
703 table3gpp->m_uXpr = 8;
704 table3gpp->m_sigXpr = 3;
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_UMi_NLOS[row][column];
712 }
713 }
714 }
715 else //(o2i)
716 {
717 table3gpp->m_numOfCluster = 12;
718 table3gpp->m_raysPerCluster = 20;
719 table3gpp->m_uLgDS = -6.62;
720 table3gpp->m_sigLgDS = 0.32;
721 table3gpp->m_uLgASD = 1.25;
722 table3gpp->m_sigLgASD = 0.42;
723 table3gpp->m_uLgASA = 1.76;
724 table3gpp->m_sigLgASA = 0.16;
725 table3gpp->m_uLgZSA = 1.01;
726 table3gpp->m_sigLgZSA = 0.43;
727 table3gpp->m_uLgZSD = uLgZSD;
728 table3gpp->m_sigLgZSD = 0.35;
729 table3gpp->m_offsetZOD = offsetZOD;
730 table3gpp->m_cDS = 11e-9;
731 table3gpp->m_cASD = 5;
732 table3gpp->m_cASA = 8;
733 table3gpp->m_cZSA = 3;
734 table3gpp->m_uK = 0;
735 table3gpp->m_sigK = 0;
736 table3gpp->m_rTau = 2.2;
737 table3gpp->m_uXpr = 9;
738 table3gpp->m_sigXpr = 5;
739 table3gpp->m_perClusterShadowingStd = 4;
740
741 for (uint8_t row = 0; row < 6; row++)
742 {
743 for (uint8_t column = 0; column < 6; column++)
744 {
745 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_O2I[row][column];
746 }
747 }
748 }
749 }
750 }
751 else if (m_scenario == "InH-OfficeMixed" || m_scenario == "InH-OfficeOpen")
752 {
753 NS_ASSERT_MSG(!o2i, "The indoor scenario does out support outdoor to indoor");
754 if (los)
755 {
756 table3gpp->m_numOfCluster = 15;
757 table3gpp->m_raysPerCluster = 20;
758 table3gpp->m_uLgDS = -0.01 * log10(1 + fcGHz) - 7.692;
759 table3gpp->m_sigLgDS = 0.18;
760 table3gpp->m_uLgASD = 1.60;
761 table3gpp->m_sigLgASD = 0.18;
762 table3gpp->m_uLgASA = -0.19 * log10(1 + fcGHz) + 1.781;
763 table3gpp->m_sigLgASA = 0.12 * log10(1 + fcGHz) + 0.119;
764 table3gpp->m_uLgZSA = -0.26 * log10(1 + fcGHz) + 1.44;
765 table3gpp->m_sigLgZSA = -0.04 * log10(1 + fcGHz) + 0.264;
766 table3gpp->m_uLgZSD = -1.43 * log10(1 + fcGHz) + 2.228;
767 table3gpp->m_sigLgZSD = 0.13 * log10(1 + fcGHz) + 0.30;
768 table3gpp->m_offsetZOD = 0;
769 table3gpp->m_cDS = 3.91e-9;
770 table3gpp->m_cASD = 5;
771 table3gpp->m_cASA = 8;
772 table3gpp->m_cZSA = 9;
773 table3gpp->m_uK = 7;
774 table3gpp->m_sigK = 4;
775 table3gpp->m_rTau = 3.6;
776 table3gpp->m_uXpr = 11;
777 table3gpp->m_sigXpr = 4;
778 table3gpp->m_perClusterShadowingStd = 6;
779
780 for (uint8_t row = 0; row < 7; row++)
781 {
782 for (uint8_t column = 0; column < 7; column++)
783 {
784 table3gpp->m_sqrtC[row][column] = sqrtC_office_LOS[row][column];
785 }
786 }
787 }
788 else
789 {
790 table3gpp->m_numOfCluster = 19;
791 table3gpp->m_raysPerCluster = 20;
792 table3gpp->m_uLgDS = -0.28 * log10(1 + fcGHz) - 7.173;
793 table3gpp->m_sigLgDS = 0.1 * log10(1 + fcGHz) + 0.055;
794 table3gpp->m_uLgASD = 1.62;
795 table3gpp->m_sigLgASD = 0.25;
796 table3gpp->m_uLgASA = -0.11 * log10(1 + fcGHz) + 1.863;
797 table3gpp->m_sigLgASA = 0.12 * log10(1 + fcGHz) + 0.059;
798 table3gpp->m_uLgZSA = -0.15 * log10(1 + fcGHz) + 1.387;
799 table3gpp->m_sigLgZSA = -0.09 * log10(1 + fcGHz) + 0.746;
800 table3gpp->m_uLgZSD = 1.08;
801 table3gpp->m_sigLgZSD = 0.36;
802 table3gpp->m_offsetZOD = 0;
803 table3gpp->m_cDS = 3.91e-9;
804 table3gpp->m_cASD = 5;
805 table3gpp->m_cASA = 11;
806 table3gpp->m_cZSA = 9;
807 table3gpp->m_uK = 0;
808 table3gpp->m_sigK = 0;
809 table3gpp->m_rTau = 3;
810 table3gpp->m_uXpr = 10;
811 table3gpp->m_sigXpr = 4;
812 table3gpp->m_perClusterShadowingStd = 3;
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_office_NLOS[row][column];
819 }
820 }
821 }
822 }
823 else if (m_scenario == "V2V-Urban")
824 {
825 if (channelCondition->IsLos())
826 {
827 // 3GPP mentioned that 3.91 ns should be used when the Cluster DS (cDS)
828 // entry is N/A.
829 table3gpp->m_numOfCluster = 12;
830 table3gpp->m_raysPerCluster = 20;
831 table3gpp->m_uLgDS = -0.2 * log10(1 + fcGHz) - 7.5;
832 table3gpp->m_sigLgDS = 0.1;
833 table3gpp->m_uLgASD = -0.1 * log10(1 + fcGHz) + 1.6;
834 table3gpp->m_sigLgASD = 0.1;
835 table3gpp->m_uLgASA = -0.1 * log10(1 + fcGHz) + 1.6;
836 table3gpp->m_sigLgASA = 0.1;
837 table3gpp->m_uLgZSA = -0.1 * log10(1 + fcGHz) + 0.73;
838 table3gpp->m_sigLgZSA = -0.04 * log10(1 + fcGHz) + 0.34;
839 table3gpp->m_uLgZSD = -0.1 * log10(1 + fcGHz) + 0.73;
840 table3gpp->m_sigLgZSD = -0.04 * log10(1 + fcGHz) + 0.34;
841 table3gpp->m_offsetZOD = 0;
842 table3gpp->m_cDS = 5;
843 table3gpp->m_cASD = 17;
844 table3gpp->m_cASA = 17;
845 table3gpp->m_cZSA = 7;
846 table3gpp->m_uK = 3.48;
847 table3gpp->m_sigK = 2;
848 table3gpp->m_rTau = 3;
849 table3gpp->m_uXpr = 9;
850 table3gpp->m_sigXpr = 3;
851 table3gpp->m_perClusterShadowingStd = 4;
852
853 for (uint8_t row = 0; row < 7; row++)
854 {
855 for (uint8_t column = 0; column < 7; column++)
856 {
857 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
858 }
859 }
860 }
861 else if (channelCondition->IsNlos())
862 {
863 table3gpp->m_numOfCluster = 19;
864 table3gpp->m_raysPerCluster = 20;
865 table3gpp->m_uLgDS = -0.3 * log10(1 + fcGHz) - 7;
866 table3gpp->m_sigLgDS = 0.28;
867 table3gpp->m_uLgASD = -0.08 * log10(1 + fcGHz) + 1.81;
868 table3gpp->m_sigLgASD = 0.05 * log10(1 + fcGHz) + 0.3;
869 table3gpp->m_uLgASA = -0.08 * log10(1 + fcGHz) + 1.81;
870 table3gpp->m_sigLgASA = 0.05 * log10(1 + fcGHz) + 0.3;
871 table3gpp->m_uLgZSA = -0.04 * log10(1 + fcGHz) + 0.92;
872 table3gpp->m_sigLgZSA = -0.07 * log10(1 + fcGHz) + 0.41;
873 table3gpp->m_uLgZSD = -0.04 * log10(1 + fcGHz) + 0.92;
874 table3gpp->m_sigLgZSD = -0.07 * log10(1 + fcGHz) + 0.41;
875 table3gpp->m_offsetZOD = 0;
876 table3gpp->m_cDS = 11;
877 table3gpp->m_cASD = 22;
878 table3gpp->m_cASA = 22;
879 table3gpp->m_cZSA = 7;
880 table3gpp->m_uK = 0; // N/A
881 table3gpp->m_sigK = 0; // N/A
882 table3gpp->m_rTau = 2.1;
883 table3gpp->m_uXpr = 8;
884 table3gpp->m_sigXpr = 3;
885 table3gpp->m_perClusterShadowingStd = 4;
886
887 for (uint8_t row = 0; row < 6; row++)
888 {
889 for (uint8_t column = 0; column < 6; column++)
890 {
891 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_NLOS[row][column];
892 }
893 }
894 }
895 else if (channelCondition->IsNlosv())
896 {
897 table3gpp->m_numOfCluster = 19;
898 table3gpp->m_raysPerCluster = 20;
899 table3gpp->m_uLgDS = -0.4 * log10(1 + fcGHz) - 7;
900 table3gpp->m_sigLgDS = 0.1;
901 table3gpp->m_uLgASD = -0.1 * log10(1 + fcGHz) + 1.7;
902 table3gpp->m_sigLgASD = 0.1;
903 table3gpp->m_uLgASA = -0.1 * log10(1 + fcGHz) + 1.7;
904 table3gpp->m_sigLgASA = 0.1;
905 table3gpp->m_uLgZSA = -0.04 * log10(1 + fcGHz) + 0.92;
906 table3gpp->m_sigLgZSA = -0.07 * log10(1 + fcGHz) + 0.41;
907 table3gpp->m_uLgZSD = -0.04 * log10(1 + fcGHz) + 0.92;
908 table3gpp->m_sigLgZSD = -0.07 * log10(1 + fcGHz) + 0.41;
909 table3gpp->m_offsetZOD = 0;
910 table3gpp->m_cDS = 11;
911 table3gpp->m_cASD = 22;
912 table3gpp->m_cASA = 22;
913 table3gpp->m_cZSA = 7;
914 table3gpp->m_uK = 0;
915 table3gpp->m_sigK = 4.5;
916 table3gpp->m_rTau = 2.1;
917 table3gpp->m_uXpr = 8;
918 table3gpp->m_sigXpr = 3;
919 table3gpp->m_perClusterShadowingStd = 4;
920
921 for (uint8_t row = 0; row < 6; row++)
922 {
923 for (uint8_t column = 0; column < 6; column++)
924 {
925 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
926 }
927 }
928 }
929 else
930 {
931 NS_FATAL_ERROR("Unknown channel condition");
932 }
933 }
934 else if (m_scenario == "V2V-Highway")
935 {
936 if (channelCondition->IsLos())
937 {
938 table3gpp->m_numOfCluster = 12;
939 table3gpp->m_raysPerCluster = 20;
940 table3gpp->m_uLgDS = -8.3;
941 table3gpp->m_sigLgDS = 0.2;
942 table3gpp->m_uLgASD = 1.4;
943 table3gpp->m_sigLgASD = 0.1;
944 table3gpp->m_uLgASA = 1.4;
945 table3gpp->m_sigLgASA = 0.1;
946 table3gpp->m_uLgZSA = -0.1 * log10(1 + fcGHz) + 0.73;
947 table3gpp->m_sigLgZSA = -0.04 * log10(1 + fcGHz) + 0.34;
948 table3gpp->m_uLgZSD = -0.1 * log10(1 + fcGHz) + 0.73;
949 table3gpp->m_sigLgZSD = -0.04 * log10(1 + fcGHz) + 0.34;
950 table3gpp->m_offsetZOD = 0;
951 table3gpp->m_cDS = 5;
952 table3gpp->m_cASD = 17;
953 table3gpp->m_cASA = 17;
954 table3gpp->m_cZSA = 7;
955 table3gpp->m_uK = 9;
956 table3gpp->m_sigK = 3.5;
957 table3gpp->m_rTau = 3;
958 table3gpp->m_uXpr = 9;
959 table3gpp->m_sigXpr = 3;
960 table3gpp->m_perClusterShadowingStd = 4;
961
962 for (uint8_t row = 0; row < 7; row++)
963 {
964 for (uint8_t column = 0; column < 7; column++)
965 {
966 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
967 }
968 }
969 }
970 else if (channelCondition->IsNlosv())
971 {
972 table3gpp->m_numOfCluster = 19;
973 table3gpp->m_raysPerCluster = 20;
974 table3gpp->m_uLgDS = -8.3;
975 table3gpp->m_sigLgDS = 0.3;
976 table3gpp->m_uLgASD = 1.5;
977 table3gpp->m_sigLgASD = 0.1;
978 table3gpp->m_uLgASA = 1.5;
979 table3gpp->m_sigLgASA = 0.1;
980 table3gpp->m_uLgZSA = -0.04 * log10(1 + fcGHz) + 0.92;
981 table3gpp->m_sigLgZSA = -0.07 * log10(1 + fcGHz) + 0.41;
982 table3gpp->m_uLgZSD = -0.04 * log10(1 + fcGHz) + 0.92;
983 table3gpp->m_sigLgZSD = -0.07 * log10(1 + fcGHz) + 0.41;
984 table3gpp->m_offsetZOD = 0;
985 table3gpp->m_cDS = 11;
986 table3gpp->m_cASD = 22;
987 table3gpp->m_cASA = 22;
988 table3gpp->m_cZSA = 7;
989 table3gpp->m_uK = 0;
990 table3gpp->m_sigK = 4.5;
991 table3gpp->m_rTau = 2.1;
992 table3gpp->m_uXpr = 8.0;
993 table3gpp->m_sigXpr = 3;
994 table3gpp->m_perClusterShadowingStd = 4;
995
996 for (uint8_t row = 0; row < 6; row++)
997 {
998 for (uint8_t column = 0; column < 6; column++)
999 {
1000 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_LOS[row][column];
1001 }
1002 }
1003 }
1004 else if (channelCondition->IsNlos())
1005 {
1007 "The fast fading parameters for the NLOS condition in the Highway scenario are not "
1008 "defined in TR 37.885, use the ones defined in TDoc R1-1803671 instead");
1009
1010 table3gpp->m_numOfCluster = 19;
1011 table3gpp->m_raysPerCluster = 20;
1012 table3gpp->m_uLgDS = -0.3 * log10(1 + fcGHz) - 7;
1013 table3gpp->m_sigLgDS = 0.28;
1014 table3gpp->m_uLgASD = -0.08 * log10(1 + fcGHz) + 1.81;
1015 table3gpp->m_sigLgASD = 0.05 * log10(1 + fcGHz) + 0.3;
1016 table3gpp->m_uLgASA = -0.08 * log10(1 + fcGHz) + 1.81;
1017 table3gpp->m_sigLgASA = 0.05 * log10(1 + fcGHz) + 0.3;
1018 table3gpp->m_uLgZSA = -0.04 * log10(1 + fcGHz) + 0.92;
1019 table3gpp->m_sigLgZSA = -0.07 * log10(1 + fcGHz) + 0.41;
1020 table3gpp->m_uLgZSD = -0.04 * log10(1 + fcGHz) + 0.92;
1021 table3gpp->m_sigLgZSD = -0.07 * log10(1 + fcGHz) + 0.41;
1022 table3gpp->m_offsetZOD = 0;
1023 table3gpp->m_cDS = 11;
1024 table3gpp->m_cASD = 22;
1025 table3gpp->m_cASA = 22;
1026 table3gpp->m_cZSA = 7;
1027 table3gpp->m_uK = 0; // N/A
1028 table3gpp->m_sigK = 0; // N/A
1029 table3gpp->m_rTau = 2.1;
1030 table3gpp->m_uXpr = 8;
1031 table3gpp->m_sigXpr = 3;
1032 table3gpp->m_perClusterShadowingStd = 4;
1033
1034 for (uint8_t row = 0; row < 6; row++)
1035 {
1036 for (uint8_t column = 0; column < 6; column++)
1037 {
1038 table3gpp->m_sqrtC[row][column] = sqrtC_UMi_NLOS[row][column];
1039 }
1040 }
1041 }
1042 else
1043 {
1044 NS_FATAL_ERROR("Unknown channel condition");
1045 }
1046 }
1047 else
1048 {
1049 NS_FATAL_ERROR("unknown scenarios");
1050 }
1051
1052 return table3gpp;
1053}
1054
1055bool
1057 Ptr<const ChannelCondition> channelCondition) const
1058{
1059 NS_LOG_FUNCTION(this);
1060
1061 bool update = false;
1062
1063 // if the channel condition is different the channel has to be updated
1064 if (!channelCondition->IsEqual(channelParams->m_losCondition, channelParams->m_o2iCondition))
1065 {
1066 NS_LOG_DEBUG("Update the channel condition");
1067 update = true;
1068 }
1069
1070 // if the coherence time is over the channel has to be updated
1071 if (!m_updatePeriod.IsZero() &&
1072 Simulator::Now() - channelParams->m_generatedTime > m_updatePeriod)
1073 {
1074 NS_LOG_DEBUG("Generation time " << channelParams->m_generatedTime.As(Time::NS) << " now "
1075 << Now().As(Time::NS));
1076 update = true;
1077 }
1078
1079 return update;
1080}
1081
1082bool
1084 Ptr<const ChannelMatrix> channelMatrix)
1085{
1086 if (channelParams->m_generatedTime > channelMatrix->m_generatedTime)
1087 {
1088 return true;
1089 }
1090 else
1091 {
1092 return false;
1093 }
1094}
1095
1101{
1102 NS_LOG_FUNCTION(this);
1103
1104 // Compute the channel params key. The key is reciprocal, i.e., key (a, b) = key (b, a)
1105 uint64_t channelParamsKey =
1106 GetKey(aMob->GetObject<Node>()->GetId(), bMob->GetObject<Node>()->GetId());
1107 // Compute the channel matrix key. The key is reciprocal, i.e., key (a, b) = key (b, a)
1108 uint64_t channelMatrixKey = GetKey(aAntenna->GetId(), bAntenna->GetId());
1109
1110 // retrieve the channel condition
1111 Ptr<const ChannelCondition> condition =
1112 m_channelConditionModel->GetChannelCondition(aMob, bMob);
1113
1114 // Check if the channel is present in the map and return it, otherwise
1115 // generate a new channel
1116 bool updateParams = false;
1117 bool updateMatrix = false;
1118 bool notFoundParams = false;
1119 bool notFoundMatrix = false;
1120 Ptr<ChannelMatrix> channelMatrix;
1121 Ptr<ThreeGppChannelParams> channelParams;
1122
1123 if (m_channelParamsMap.find(channelParamsKey) != m_channelParamsMap.end())
1124 {
1125 channelParams = m_channelParamsMap[channelParamsKey];
1126 // check if it has to be updated
1127 updateParams = ChannelParamsNeedsUpdate(channelParams, condition);
1128 }
1129 else
1130 {
1131 NS_LOG_DEBUG("channel params not found");
1132 notFoundParams = true;
1133 }
1134
1135 double x = aMob->GetPosition().x - bMob->GetPosition().x;
1136 double y = aMob->GetPosition().y - bMob->GetPosition().y;
1137 double distance2D = sqrt(x * x + y * y);
1138
1139 // NOTE we assume hUT = min (height(a), height(b)) and
1140 // hBS = max (height (a), height (b))
1141 double hUt = std::min(aMob->GetPosition().z, bMob->GetPosition().z);
1142 double hBs = std::max(aMob->GetPosition().z, bMob->GetPosition().z);
1143
1144 // get the 3GPP parameters
1145 Ptr<const ParamsTable> table3gpp = GetThreeGppTable(condition, hBs, hUt, distance2D);
1146
1147 if (notFoundParams || updateParams)
1148 {
1149 // Step 4: Generate large scale parameters. All LSPS are uncorrelated.
1150 // Step 5: Generate Delays.
1151 // Step 6: Generate cluster powers.
1152 // Step 7: Generate arrival and departure angles for both azimuth and elevation.
1153 // Step 8: Coupling of rays within a cluster for both azimuth and elevation
1154 // shuffle all the arrays to perform random coupling
1155 // Step 9: Generate the cross polarization power ratios
1156 // Step 10: Draw initial phases
1157 channelParams = GenerateChannelParameters(condition, table3gpp, aMob, bMob);
1158 // store or replace the channel parameters
1159 m_channelParamsMap[channelParamsKey] = channelParams;
1160 }
1161
1162 if (m_channelMatrixMap.find(channelMatrixKey) != m_channelMatrixMap.end())
1163 {
1164 // channel matrix present in the map
1165 NS_LOG_DEBUG("channel matrix present in the map");
1166 channelMatrix = m_channelMatrixMap[channelMatrixKey];
1167 updateMatrix = ChannelMatrixNeedsUpdate(channelParams, channelMatrix);
1168 }
1169 else
1170 {
1171 NS_LOG_DEBUG("channel matrix not found");
1172 notFoundMatrix = true;
1173 }
1174
1175 // If the channel is not present in the map or if it has to be updated
1176 // generate a new realization
1177 if (notFoundMatrix || updateMatrix)
1178 {
1179 // channel matrix not found or has to be updated, generate a new one
1180 channelMatrix = GetNewChannel(channelParams, table3gpp, aMob, bMob, aAntenna, bAntenna);
1181 channelMatrix->m_antennaPair =
1182 std::make_pair(aAntenna->GetId(),
1183 bAntenna->GetId()); // save antenna pair, with the exact order of s and u
1184 // antennas at the moment of the channel generation
1185
1186 // store or replace the channel matrix in the channel map
1187 m_channelMatrixMap[channelMatrixKey] = channelMatrix;
1188 }
1189
1190 return channelMatrix;
1191}
1192
1195{
1196 NS_LOG_FUNCTION(this);
1197
1198 // Compute the channel key. The key is reciprocal, i.e., key (a, b) = key (b, a)
1199 uint64_t channelParamsKey =
1200 GetKey(aMob->GetObject<Node>()->GetId(), bMob->GetObject<Node>()->GetId());
1201
1202 if (m_channelParamsMap.find(channelParamsKey) != m_channelParamsMap.end())
1203 {
1204 return m_channelParamsMap.find(channelParamsKey)->second;
1205 }
1206 else
1207 {
1208 NS_LOG_WARN("Channel params map not found. Returning a nullptr.");
1209 return nullptr;
1210 }
1211}
1212
1215 const Ptr<const ParamsTable> table3gpp,
1216 const Ptr<const MobilityModel> aMob,
1217 const Ptr<const MobilityModel> bMob) const
1218{
1219 NS_LOG_FUNCTION(this);
1220 // create a channel matrix instance
1221 Ptr<ThreeGppChannelParams> channelParams = Create<ThreeGppChannelParams>();
1222 channelParams->m_generatedTime = Simulator::Now();
1223 channelParams->m_nodeIds =
1224 std::make_pair(aMob->GetObject<Node>()->GetId(), bMob->GetObject<Node>()->GetId());
1225 channelParams->m_losCondition = channelCondition->GetLosCondition();
1226 channelParams->m_o2iCondition = channelCondition->GetO2iCondition();
1227
1228 // Step 4: Generate large scale parameters. All LSPS are uncorrelated.
1229 DoubleVector LSPsIndep;
1230 DoubleVector LSPs;
1231 uint8_t paramNum = 6;
1232 if (channelParams->m_losCondition == ChannelCondition::LOS)
1233 {
1234 paramNum = 7;
1235 }
1236
1237 // Generate paramNum independent LSPs.
1238 for (uint8_t iter = 0; iter < paramNum; iter++)
1239 {
1240 LSPsIndep.push_back(m_normalRv->GetValue());
1241 }
1242 for (uint8_t row = 0; row < paramNum; row++)
1243 {
1244 double temp = 0;
1245 for (uint8_t column = 0; column < paramNum; column++)
1246 {
1247 temp += table3gpp->m_sqrtC[row][column] * LSPsIndep[column];
1248 }
1249 LSPs.push_back(temp);
1250 }
1251
1252 // NOTE the shadowing is generated in the propagation loss model
1253 double DS;
1254 double ASD;
1255 double ASA;
1256 double ZSA;
1257 double ZSD;
1258 double kFactor = 0;
1259 if (channelParams->m_losCondition == ChannelCondition::LOS)
1260 {
1261 kFactor = LSPs[1] * table3gpp->m_sigK + table3gpp->m_uK;
1262 DS = pow(10, LSPs[2] * table3gpp->m_sigLgDS + table3gpp->m_uLgDS);
1263 ASD = pow(10, LSPs[3] * table3gpp->m_sigLgASD + table3gpp->m_uLgASD);
1264 ASA = pow(10, LSPs[4] * table3gpp->m_sigLgASA + table3gpp->m_uLgASA);
1265 ZSD = pow(10, LSPs[5] * table3gpp->m_sigLgZSD + table3gpp->m_uLgZSD);
1266 ZSA = pow(10, LSPs[6] * table3gpp->m_sigLgZSA + table3gpp->m_uLgZSA);
1267 }
1268 else
1269 {
1270 DS = pow(10, LSPs[1] * table3gpp->m_sigLgDS + table3gpp->m_uLgDS);
1271 ASD = pow(10, LSPs[2] * table3gpp->m_sigLgASD + table3gpp->m_uLgASD);
1272 ASA = pow(10, LSPs[3] * table3gpp->m_sigLgASA + table3gpp->m_uLgASA);
1273 ZSD = pow(10, LSPs[4] * table3gpp->m_sigLgZSD + table3gpp->m_uLgZSD);
1274 ZSA = pow(10, LSPs[5] * table3gpp->m_sigLgZSA + table3gpp->m_uLgZSA);
1275 }
1276 ASD = std::min(ASD, 104.0);
1277 ASA = std::min(ASA, 104.0);
1278 ZSD = std::min(ZSD, 52.0);
1279 ZSA = std::min(ZSA, 52.0);
1280
1281 // save DS and K_factor parameters in the structure
1282 channelParams->m_DS = DS;
1283 channelParams->m_K_factor = kFactor;
1284
1285 NS_LOG_INFO("K-factor=" << kFactor << ", DS=" << DS << ", ASD=" << ASD << ", ASA=" << ASA
1286 << ", ZSD=" << ZSD << ", ZSA=" << ZSA);
1287
1288 // Step 5: Generate Delays.
1289 DoubleVector clusterDelay;
1290 double minTau = 100.0;
1291 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1292 {
1293 double tau = -1 * table3gpp->m_rTau * DS * log(m_uniformRv->GetValue(0, 1)); //(7.5-1)
1294 if (minTau > tau)
1295 {
1296 minTau = tau;
1297 }
1298 clusterDelay.push_back(tau);
1299 }
1300
1301 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1302 {
1303 clusterDelay[cIndex] -= minTau;
1304 }
1305 std::sort(clusterDelay.begin(), clusterDelay.end()); //(7.5-2)
1306
1307 /* since the scaled Los delays are not to be used in cluster power generation,
1308 * we will generate cluster power first and resume to compute Los cluster delay later.*/
1309
1310 // Step 6: Generate cluster powers.
1311 DoubleVector clusterPower;
1312 double powerSum = 0;
1313 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1314 {
1315 double power =
1316 exp(-1 * clusterDelay[cIndex] * (table3gpp->m_rTau - 1) / table3gpp->m_rTau / DS) *
1317 pow(10,
1318 -1 * m_normalRv->GetValue() * table3gpp->m_perClusterShadowingStd / 10.0); //(7.5-5)
1319 powerSum += power;
1320 clusterPower.push_back(power);
1321 }
1322 channelParams->m_clusterPower = clusterPower;
1323
1324 double powerMax = 0;
1325
1326 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1327 {
1328 channelParams->m_clusterPower[cIndex] =
1329 channelParams->m_clusterPower[cIndex] / powerSum; //(7.5-6)
1330 }
1331
1332 DoubleVector clusterPowerForAngles; // this power is only for equation (7.5-9) and (7.5-14), not
1333 // for (7.5-22)
1334 if (channelParams->m_losCondition == ChannelCondition::LOS)
1335 {
1336 double kLinear = pow(10, kFactor / 10.0);
1337
1338 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1339 {
1340 if (cIndex == 0)
1341 {
1342 clusterPowerForAngles.push_back(channelParams->m_clusterPower[cIndex] /
1343 (1 + kLinear) +
1344 kLinear / (1 + kLinear)); //(7.5-8)
1345 }
1346 else
1347 {
1348 clusterPowerForAngles.push_back(channelParams->m_clusterPower[cIndex] /
1349 (1 + kLinear)); //(7.5-8)
1350 }
1351 if (powerMax < clusterPowerForAngles[cIndex])
1352 {
1353 powerMax = clusterPowerForAngles[cIndex];
1354 }
1355 }
1356 }
1357 else
1358 {
1359 for (uint8_t cIndex = 0; cIndex < table3gpp->m_numOfCluster; cIndex++)
1360 {
1361 clusterPowerForAngles.push_back(channelParams->m_clusterPower[cIndex]); //(7.5-6)
1362 if (powerMax < clusterPowerForAngles[cIndex])
1363 {
1364 powerMax = clusterPowerForAngles[cIndex];
1365 }
1366 }
1367 }
1368
1369 // remove clusters with less than -25 dB power compared to the maxim cluster power;
1370 // double thresh = pow(10, -2.5);
1371 double thresh = 0.0032;
1372 for (uint8_t cIndex = table3gpp->m_numOfCluster; cIndex > 0; cIndex--)
1373 {
1374 if (clusterPowerForAngles[cIndex - 1] < thresh * powerMax)
1375 {
1376 clusterPowerForAngles.erase(clusterPowerForAngles.begin() + cIndex - 1);
1377 channelParams->m_clusterPower.erase(channelParams->m_clusterPower.begin() + cIndex - 1);
1378 clusterDelay.erase(clusterDelay.begin() + cIndex - 1);
1379 }
1380 }
1381
1382 NS_ASSERT(channelParams->m_clusterPower.size() < UINT8_MAX);
1383 channelParams->m_reducedClusterNumber = channelParams->m_clusterPower.size();
1384 // Resume step 5 to compute the delay for LoS condition.
1385 if (channelParams->m_losCondition == ChannelCondition::LOS)
1386 {
1387 double cTau =
1388 0.7705 - 0.0433 * kFactor + 2e-4 * pow(kFactor, 2) + 17e-6 * pow(kFactor, 3); //(7.5-3)
1389 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1390 {
1391 clusterDelay[cIndex] = clusterDelay[cIndex] / cTau; //(7.5-4)
1392 }
1393 }
1394
1395 // Step 7: Generate arrival and departure angles for both azimuth and elevation.
1396
1397 double cNlos;
1398 // According to table 7.5-6, only cluster number equals to 8, 10, 11, 12, 19 and 20 is valid.
1399 // Not sure why the other cases are in Table 7.5-2.
1400 switch (table3gpp->m_numOfCluster) // Table 7.5-2
1401 {
1402 case 4:
1403 cNlos = 0.779;
1404 break;
1405 case 5:
1406 cNlos = 0.860;
1407 break;
1408 case 8:
1409 cNlos = 1.018;
1410 break;
1411 case 10:
1412 cNlos = 1.090;
1413 break;
1414 case 11:
1415 cNlos = 1.123;
1416 break;
1417 case 12:
1418 cNlos = 1.146;
1419 break;
1420 case 14:
1421 cNlos = 1.190;
1422 break;
1423 case 15:
1424 cNlos = 1.221;
1425 break;
1426 case 16:
1427 cNlos = 1.226;
1428 break;
1429 case 19:
1430 cNlos = 1.273;
1431 break;
1432 case 20:
1433 cNlos = 1.289;
1434 break;
1435 default:
1436 NS_FATAL_ERROR("Invalid cluster number");
1437 }
1438
1439 double cPhi = cNlos;
1440
1441 if (channelParams->m_losCondition == ChannelCondition::LOS)
1442 {
1443 cPhi *= (1.1035 - 0.028 * kFactor - 2e-3 * pow(kFactor, 2) +
1444 1e-4 * pow(kFactor, 3)); //(7.5-10))
1445 }
1446
1447 switch (table3gpp->m_numOfCluster) // Table 7.5-4
1448 {
1449 case 8:
1450 cNlos = 0.889;
1451 break;
1452 case 10:
1453 cNlos = 0.957;
1454 break;
1455 case 11:
1456 cNlos = 1.031;
1457 break;
1458 case 12:
1459 cNlos = 1.104;
1460 break;
1461 case 15:
1462 cNlos = 1.1088;
1463 break;
1464 case 19:
1465 cNlos = 1.184;
1466 break;
1467 case 20:
1468 cNlos = 1.178;
1469 break;
1470 default:
1471 NS_FATAL_ERROR("Invalid cluster number");
1472 }
1473
1474 double cTheta = cNlos;
1475 if (channelCondition->IsLos())
1476 {
1477 cTheta *= (1.3086 + 0.0339 * kFactor - 0.0077 * pow(kFactor, 2) +
1478 2e-4 * pow(kFactor, 3)); //(7.5-15)
1479 }
1480
1481 DoubleVector clusterAoa;
1482 DoubleVector clusterAod;
1483 DoubleVector clusterZoa;
1484 DoubleVector clusterZod;
1485 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1486 {
1487 double logCalc = -1 * log(clusterPowerForAngles[cIndex] / powerMax);
1488 double angle = 2 * sqrt(logCalc) / 1.4 / cPhi; //(7.5-9)
1489 clusterAoa.push_back(ASA * angle);
1490 clusterAod.push_back(ASD * angle);
1491 angle = logCalc / cTheta; //(7.5-14)
1492 clusterZoa.push_back(ZSA * angle);
1493 clusterZod.push_back(ZSD * angle);
1494 }
1495
1496 Angles sAngle(bMob->GetPosition(), aMob->GetPosition());
1497 Angles uAngle(aMob->GetPosition(), bMob->GetPosition());
1498
1499 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1500 {
1501 int Xn = 1;
1502 if (m_uniformRv->GetValue(0, 1) < 0.5)
1503 {
1504 Xn = -1;
1505 }
1506 clusterAoa[cIndex] = clusterAoa[cIndex] * Xn + (m_normalRv->GetValue() * ASA / 7.0) +
1507 RadiansToDegrees(uAngle.GetAzimuth()); //(7.5-11)
1508 clusterAod[cIndex] = clusterAod[cIndex] * Xn + (m_normalRv->GetValue() * ASD / 7.0) +
1509 RadiansToDegrees(sAngle.GetAzimuth());
1510 if (channelCondition->IsO2i())
1511 {
1512 clusterZoa[cIndex] =
1513 clusterZoa[cIndex] * Xn + (m_normalRv->GetValue() * ZSA / 7.0) + 90; //(7.5-16)
1514 }
1515 else
1516 {
1517 clusterZoa[cIndex] = clusterZoa[cIndex] * Xn + (m_normalRv->GetValue() * ZSA / 7.0) +
1518 RadiansToDegrees(uAngle.GetInclination()); //(7.5-16)
1519 }
1520 clusterZod[cIndex] = clusterZod[cIndex] * Xn + (m_normalRv->GetValue() * ZSD / 7.0) +
1522 table3gpp->m_offsetZOD; //(7.5-19)
1523 }
1524
1525 if (channelParams->m_losCondition == ChannelCondition::LOS)
1526 {
1527 // The 7.5-12 can be rewrite as Theta_n,ZOA = Theta_n,ZOA - (Theta_1,ZOA - Theta_LOS,ZOA) =
1528 // Theta_n,ZOA - diffZOA, Similar as AOD, ZSA and ZSD.
1529 double diffAoa = clusterAoa[0] - RadiansToDegrees(uAngle.GetAzimuth());
1530 double diffAod = clusterAod[0] - RadiansToDegrees(sAngle.GetAzimuth());
1531 double diffZsa = clusterZoa[0] - RadiansToDegrees(uAngle.GetInclination());
1532 double diffZsd = clusterZod[0] - RadiansToDegrees(sAngle.GetInclination());
1533
1534 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1535 {
1536 clusterAoa[cIndex] -= diffAoa; //(7.5-12)
1537 clusterAod[cIndex] -= diffAod;
1538 clusterZoa[cIndex] -= diffZsa; //(7.5-17)
1539 clusterZod[cIndex] -= diffZsd;
1540 }
1541 }
1542
1543 double sizeTemp = clusterZoa.size();
1544 for (uint8_t ind = 0; ind < 4; ind++)
1545 {
1546 DoubleVector angleDegree;
1547 switch (ind)
1548 {
1549 case 0:
1550 angleDegree = clusterAoa;
1551 break;
1552 case 1:
1553 angleDegree = clusterZoa;
1554 break;
1555 case 2:
1556 angleDegree = clusterAod;
1557 break;
1558 case 3:
1559 angleDegree = clusterZod;
1560 break;
1561 default:
1562 NS_FATAL_ERROR("Programming Error");
1563 }
1564 for (uint8_t nIndex = 0; nIndex < sizeTemp; nIndex++)
1565 {
1566 while (angleDegree[nIndex] > 360)
1567 {
1568 angleDegree[nIndex] -= 360;
1569 }
1570
1571 while (angleDegree[nIndex] < 0)
1572 {
1573 angleDegree[nIndex] += 360;
1574 }
1575
1576 if (ind == 1 || ind == 3)
1577 {
1578 if (angleDegree[nIndex] > 180)
1579 {
1580 angleDegree[nIndex] = 360 - angleDegree[nIndex];
1581 }
1582 }
1583 }
1584 switch (ind)
1585 {
1586 case 0:
1587 clusterAoa = angleDegree;
1588 break;
1589 case 1:
1590 clusterZoa = angleDegree;
1591 break;
1592 case 2:
1593 clusterAod = angleDegree;
1594 break;
1595 case 3:
1596 clusterZod = angleDegree;
1597 break;
1598 default:
1599 NS_FATAL_ERROR("Programming Error");
1600 }
1601 }
1602
1603 DoubleVector attenuationDb;
1604 if (m_blockage)
1605 {
1606 attenuationDb = CalcAttenuationOfBlockage(channelParams, clusterAoa, clusterZoa);
1607 for (uint8_t cInd = 0; cInd < channelParams->m_reducedClusterNumber; cInd++)
1608 {
1609 channelParams->m_clusterPower[cInd] =
1610 channelParams->m_clusterPower[cInd] / pow(10, attenuationDb[cInd] / 10.0);
1611 }
1612 }
1613 else
1614 {
1615 attenuationDb.push_back(0);
1616 }
1617
1618 // store attenuation
1619 channelParams->m_attenuation_dB = attenuationDb;
1620
1621 // Step 8: Coupling of rays within a cluster for both azimuth and elevation
1622 // shuffle all the arrays to perform random coupling
1624 channelParams->m_reducedClusterNumber,
1625 DoubleVector(table3gpp->m_raysPerCluster,
1626 0)); // rayAoaRadian[n][m], where n is cluster index, m is ray index
1628 channelParams->m_reducedClusterNumber,
1629 DoubleVector(table3gpp->m_raysPerCluster,
1630 0)); // rayAodRadian[n][m], where n is cluster index, m is ray index
1632 channelParams->m_reducedClusterNumber,
1633 DoubleVector(table3gpp->m_raysPerCluster,
1634 0)); // rayZoaRadian[n][m], where n is cluster index, m is ray index
1636 channelParams->m_reducedClusterNumber,
1637 DoubleVector(table3gpp->m_raysPerCluster,
1638 0)); // rayZodRadian[n][m], where n is cluster index, m is ray index
1639
1640 for (uint8_t nInd = 0; nInd < channelParams->m_reducedClusterNumber; nInd++)
1641 {
1642 for (uint8_t mInd = 0; mInd < table3gpp->m_raysPerCluster; mInd++)
1643 {
1644 double tempAoa = clusterAoa[nInd] + table3gpp->m_cASA * offSetAlpha[mInd]; //(7.5-13)
1645 double tempZoa = clusterZoa[nInd] + table3gpp->m_cZSA * offSetAlpha[mInd]; //(7.5-18)
1646 std::tie(rayAoaRadian[nInd][mInd], rayZoaRadian[nInd][mInd]) =
1647 WrapAngles(DegreesToRadians(tempAoa), DegreesToRadians(tempZoa));
1648
1649 double tempAod = clusterAod[nInd] + table3gpp->m_cASD * offSetAlpha[mInd]; //(7.5-13)
1650 double tempZod = clusterZod[nInd] +
1651 0.375 * pow(10, table3gpp->m_uLgZSD) * offSetAlpha[mInd]; //(7.5-20)
1652 std::tie(rayAodRadian[nInd][mInd], rayZodRadian[nInd][mInd]) =
1653 WrapAngles(DegreesToRadians(tempAod), DegreesToRadians(tempZod));
1654 }
1655 }
1656
1657 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1658 {
1659 Shuffle(&rayAodRadian[cIndex][0], &rayAodRadian[cIndex][table3gpp->m_raysPerCluster]);
1660 Shuffle(&rayAoaRadian[cIndex][0], &rayAoaRadian[cIndex][table3gpp->m_raysPerCluster]);
1661 Shuffle(&rayZodRadian[cIndex][0], &rayZodRadian[cIndex][table3gpp->m_raysPerCluster]);
1662 Shuffle(&rayZoaRadian[cIndex][0], &rayZoaRadian[cIndex][table3gpp->m_raysPerCluster]);
1663 }
1664
1665 // store values
1666 channelParams->m_rayAodRadian = rayAodRadian;
1667 channelParams->m_rayAoaRadian = rayAoaRadian;
1668 channelParams->m_rayZodRadian = rayZodRadian;
1669 channelParams->m_rayZoaRadian = rayZoaRadian;
1670
1671 // Step 9: Generate the cross polarization power ratios
1672 // Step 10: Draw initial phases
1673 Double2DVector crossPolarizationPowerRatios; // vector containing the cross polarization power
1674 // ratios, as defined by 7.5-21
1675 Double3DVector clusterPhase; // rayAoaRadian[n][m], where n is cluster index, m is ray index
1676 for (uint8_t nInd = 0; nInd < channelParams->m_reducedClusterNumber; nInd++)
1677 {
1678 DoubleVector temp; // used to store the XPR values
1680 temp2; // used to store the PHI values for all the possible combination of polarization
1681 for (uint8_t mInd = 0; mInd < table3gpp->m_raysPerCluster; mInd++)
1682 {
1683 double uXprLinear = pow(10, table3gpp->m_uXpr / 10.0); // convert to linear
1684 double sigXprLinear = pow(10, table3gpp->m_sigXpr / 10.0); // convert to linear
1685
1686 temp.push_back(
1687 std::pow(10, (m_normalRv->GetValue() * sigXprLinear + uXprLinear) / 10.0));
1688 DoubleVector temp3; // used to store the PHI values
1689 for (uint8_t pInd = 0; pInd < 4; pInd++)
1690 {
1691 temp3.push_back(m_uniformRv->GetValue(-1 * M_PI, M_PI));
1692 }
1693 temp2.push_back(temp3);
1694 }
1695 crossPolarizationPowerRatios.push_back(temp);
1696 clusterPhase.push_back(temp2);
1697 }
1698 // store the cluster phase
1699 channelParams->m_clusterPhase = clusterPhase;
1700 channelParams->m_crossPolarizationPowerRatios = crossPolarizationPowerRatios;
1701
1702 uint8_t cluster1st = 0;
1703 uint8_t cluster2nd = 0; // first and second strongest cluster;
1704 double maxPower = 0;
1705 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1706 {
1707 if (maxPower < channelParams->m_clusterPower[cIndex])
1708 {
1709 maxPower = channelParams->m_clusterPower[cIndex];
1710 cluster1st = cIndex;
1711 }
1712 }
1713 channelParams->m_cluster1st = cluster1st;
1714 maxPower = 0;
1715 for (uint8_t cIndex = 0; cIndex < channelParams->m_reducedClusterNumber; cIndex++)
1716 {
1717 if (maxPower < channelParams->m_clusterPower[cIndex] && cluster1st != cIndex)
1718 {
1719 maxPower = channelParams->m_clusterPower[cIndex];
1720 cluster2nd = cIndex;
1721 }
1722 }
1723 channelParams->m_cluster2nd = cluster2nd;
1724
1725 NS_LOG_INFO("1st strongest cluster:" << +cluster1st
1726 << ", 2nd strongest cluster:" << +cluster2nd);
1727
1728 // store the delays and the angles for the subclusters
1729 if (cluster1st == cluster2nd)
1730 {
1731 clusterDelay.push_back(clusterDelay[cluster1st] + 1.28 * table3gpp->m_cDS);
1732 clusterDelay.push_back(clusterDelay[cluster1st] + 2.56 * table3gpp->m_cDS);
1733
1734 clusterAoa.push_back(clusterAoa[cluster1st]);
1735 clusterAoa.push_back(clusterAoa[cluster1st]);
1736
1737 clusterZoa.push_back(clusterZoa[cluster1st]);
1738 clusterZoa.push_back(clusterZoa[cluster1st]);
1739
1740 clusterAod.push_back(clusterAod[cluster1st]);
1741 clusterAod.push_back(clusterAod[cluster1st]);
1742
1743 clusterZod.push_back(clusterZod[cluster1st]);
1744 clusterZod.push_back(clusterZod[cluster1st]);
1745 }
1746 else
1747 {
1748 double min;
1749 double max;
1750 if (cluster1st < cluster2nd)
1751 {
1752 min = cluster1st;
1753 max = cluster2nd;
1754 }
1755 else
1756 {
1757 min = cluster2nd;
1758 max = cluster1st;
1759 }
1760 clusterDelay.push_back(clusterDelay[min] + 1.28 * table3gpp->m_cDS);
1761 clusterDelay.push_back(clusterDelay[min] + 2.56 * table3gpp->m_cDS);
1762 clusterDelay.push_back(clusterDelay[max] + 1.28 * table3gpp->m_cDS);
1763 clusterDelay.push_back(clusterDelay[max] + 2.56 * table3gpp->m_cDS);
1764
1765 clusterAoa.push_back(clusterAoa[min]);
1766 clusterAoa.push_back(clusterAoa[min]);
1767 clusterAoa.push_back(clusterAoa[max]);
1768 clusterAoa.push_back(clusterAoa[max]);
1769
1770 clusterZoa.push_back(clusterZoa[min]);
1771 clusterZoa.push_back(clusterZoa[min]);
1772 clusterZoa.push_back(clusterZoa[max]);
1773 clusterZoa.push_back(clusterZoa[max]);
1774
1775 clusterAod.push_back(clusterAod[min]);
1776 clusterAod.push_back(clusterAod[min]);
1777 clusterAod.push_back(clusterAod[max]);
1778 clusterAod.push_back(clusterAod[max]);
1779
1780 clusterZod.push_back(clusterZod[min]);
1781 clusterZod.push_back(clusterZod[min]);
1782 clusterZod.push_back(clusterZod[max]);
1783 clusterZod.push_back(clusterZod[max]);
1784 }
1785
1786 channelParams->m_delay = clusterDelay;
1787 channelParams->m_angle.clear();
1788 channelParams->m_angle.push_back(clusterAoa);
1789 channelParams->m_angle.push_back(clusterZoa);
1790 channelParams->m_angle.push_back(clusterAod);
1791 channelParams->m_angle.push_back(clusterZod);
1792
1793 // Compute alpha and D as described in 3GPP TR 37.885 v15.3.0, Sec. 6.2.3
1794 // These terms account for an additional Doppler contribution due to the
1795 // presence of moving objects in the surrounding environment, such as in
1796 // vehicular scenarios.
1797 // This contribution is applied only to the delayed (reflected) paths and
1798 // must be properly configured by setting the value of
1799 // m_vScatt, which is defined as "maximum speed of the vehicle in the
1800 // layout".
1801 // By default, m_vScatt is set to 0, so there is no additional Doppler
1802 // contribution.
1803
1804 DoubleVector dopplerTermAlpha;
1805 DoubleVector dopplerTermD;
1806
1807 // 2 or 4 is added to account for additional subrays for the 1st and 2nd clusters, if there is
1808 // only one cluster then would be added 2 more subrays (see creation of Husn channel matrix)
1809 uint8_t updatedClusterNumber = (channelParams->m_reducedClusterNumber == 1)
1810 ? channelParams->m_reducedClusterNumber + 2
1811 : channelParams->m_reducedClusterNumber + 4;
1812
1813 for (uint8_t cIndex = 0; cIndex < updatedClusterNumber; cIndex++)
1814 {
1815 double alpha = 0;
1816 double D = 0;
1817 if (cIndex != 0)
1818 {
1819 alpha = m_uniformRvDoppler->GetValue(-1, 1);
1821 }
1822 dopplerTermAlpha.push_back(alpha);
1823 dopplerTermD.push_back(D);
1824 }
1825 channelParams->m_alpha = dopplerTermAlpha;
1826 channelParams->m_D = dopplerTermD;
1827
1828 return channelParams;
1829}
1830
1833 Ptr<const ParamsTable> table3gpp,
1834 const Ptr<const MobilityModel> sMob,
1835 const Ptr<const MobilityModel> uMob,
1837 Ptr<const PhasedArrayModel> uAntenna) const
1838{
1839 NS_LOG_FUNCTION(this);
1840
1841 NS_ASSERT_MSG(m_frequency > 0.0, "Set the operating frequency first!");
1842
1843 // create a channel matrix instance
1844 Ptr<ChannelMatrix> channelMatrix = Create<ChannelMatrix>();
1845 channelMatrix->m_generatedTime = Simulator::Now();
1846 // save in which order is generated this matrix
1847 channelMatrix->m_nodeIds =
1848 std::make_pair(sMob->GetObject<Node>()->GetId(), uMob->GetObject<Node>()->GetId());
1849 // check if channelParams structure is generated in direction s-to-u or u-to-s
1850 bool isSameDirection = (channelParams->m_nodeIds == channelMatrix->m_nodeIds);
1851
1856
1857 // if channel params is generated in the same direction in which we
1858 // generate the channel matrix, angles and zenith od departure and arrival are ok,
1859 // just set them to corresponding variable that will be used for the generation
1860 // of channel matrix, otherwise we need to flip angles and zeniths of departure and arrival
1861 if (isSameDirection)
1862 {
1863 rayAodRadian = channelParams->m_rayAodRadian;
1864 rayAoaRadian = channelParams->m_rayAoaRadian;
1865 rayZodRadian = channelParams->m_rayZodRadian;
1866 rayZoaRadian = channelParams->m_rayZoaRadian;
1867 }
1868 else
1869 {
1870 rayAodRadian = channelParams->m_rayAoaRadian;
1871 rayAoaRadian = channelParams->m_rayAodRadian;
1872 rayZodRadian = channelParams->m_rayZoaRadian;
1873 rayZoaRadian = channelParams->m_rayZodRadian;
1874 }
1875
1876 // Step 11: Generate channel coefficients for each cluster n and each receiver
1877 // and transmitter element pair u,s.
1878 // where n is cluster index, u and s are receive and transmit antenna element.
1879 size_t uSize = uAntenna->GetNumberOfElements();
1880 size_t sSize = sAntenna->GetNumberOfElements();
1881
1882 // NOTE: Since each of the strongest 2 clusters are divided into 3 sub-clusters,
1883 // the total cluster will generally be numReducedCLuster + 4.
1884 // However, it might be that m_cluster1st = m_cluster2nd. In this case the
1885 // total number of clusters will be numReducedCLuster + 2.
1886 uint16_t numOverallCluster = (channelParams->m_cluster1st != channelParams->m_cluster2nd)
1887 ? channelParams->m_reducedClusterNumber + 4
1888 : channelParams->m_reducedClusterNumber + 2;
1889 Complex3DVector hUsn(uSize, sSize, numOverallCluster); // channel coefficient hUsn (u, s, n);
1890 NS_ASSERT(channelParams->m_reducedClusterNumber <= channelParams->m_clusterPhase.size());
1891 NS_ASSERT(channelParams->m_reducedClusterNumber <= channelParams->m_clusterPower.size());
1892 NS_ASSERT(channelParams->m_reducedClusterNumber <=
1893 channelParams->m_crossPolarizationPowerRatios.size());
1894 NS_ASSERT(channelParams->m_reducedClusterNumber <= rayZoaRadian.size());
1895 NS_ASSERT(channelParams->m_reducedClusterNumber <= rayZodRadian.size());
1896 NS_ASSERT(channelParams->m_reducedClusterNumber <= rayAoaRadian.size());
1897 NS_ASSERT(channelParams->m_reducedClusterNumber <= rayAodRadian.size());
1898 NS_ASSERT(table3gpp->m_raysPerCluster <= channelParams->m_clusterPhase[0].size());
1899 NS_ASSERT(table3gpp->m_raysPerCluster <=
1900 channelParams->m_crossPolarizationPowerRatios[0].size());
1901 NS_ASSERT(table3gpp->m_raysPerCluster <= rayZoaRadian[0].size());
1902 NS_ASSERT(table3gpp->m_raysPerCluster <= rayZodRadian[0].size());
1903 NS_ASSERT(table3gpp->m_raysPerCluster <= rayAoaRadian[0].size());
1904 NS_ASSERT(table3gpp->m_raysPerCluster <= rayAodRadian[0].size());
1905
1906 double x = sMob->GetPosition().x - uMob->GetPosition().x;
1907 double y = sMob->GetPosition().y - uMob->GetPosition().y;
1908 double distance2D = sqrt(x * x + y * y);
1909 // NOTE we assume hUT = min (height(a), height(b)) and
1910 // hBS = max (height (a), height (b))
1911 double hUt = std::min(sMob->GetPosition().z, uMob->GetPosition().z);
1912 double hBs = std::max(sMob->GetPosition().z, uMob->GetPosition().z);
1913 // compute the 3D distance using eq. 7.4-1
1914 double distance3D = std::sqrt(distance2D * distance2D + (hBs - hUt) * (hBs - hUt));
1915
1916 Angles sAngle(uMob->GetPosition(), sMob->GetPosition());
1917 Angles uAngle(sMob->GetPosition(), uMob->GetPosition());
1918
1919 Complex2DVector raysPreComp(channelParams->m_reducedClusterNumber,
1920 table3gpp->m_raysPerCluster); // stores part of the ray expression,
1921 // cached as independent from the u- and s-indexes
1922 Double2DVector sinCosA; // cached multiplications of sin and cos of the ZoA and AoA angles
1923 Double2DVector sinSinA; // cached multiplications of sines of the ZoA and AoA angles
1924 Double2DVector cosZoA; // cached cos of the ZoA angle
1925 Double2DVector sinCosD; // cached multiplications of sin and cos of the ZoD and AoD angles
1926 Double2DVector sinSinD; // cached multiplications of the cosines of the ZoA and AoA angles
1927 Double2DVector cosZoD; // cached cos of the ZoD angle
1928
1929 // resize to appropriate dimensions
1930 sinCosA.resize(channelParams->m_reducedClusterNumber);
1931 sinSinA.resize(channelParams->m_reducedClusterNumber);
1932 cosZoA.resize(channelParams->m_reducedClusterNumber);
1933 sinCosD.resize(channelParams->m_reducedClusterNumber);
1934 sinSinD.resize(channelParams->m_reducedClusterNumber);
1935 cosZoD.resize(channelParams->m_reducedClusterNumber);
1936 for (uint8_t nIndex = 0; nIndex < channelParams->m_reducedClusterNumber; nIndex++)
1937 {
1938 sinCosA[nIndex].resize(table3gpp->m_raysPerCluster);
1939 sinSinA[nIndex].resize(table3gpp->m_raysPerCluster);
1940 cosZoA[nIndex].resize(table3gpp->m_raysPerCluster);
1941 sinCosD[nIndex].resize(table3gpp->m_raysPerCluster);
1942 sinSinD[nIndex].resize(table3gpp->m_raysPerCluster);
1943 cosZoD[nIndex].resize(table3gpp->m_raysPerCluster);
1944 }
1945 // pre-compute the terms which are independent from uIndex and sIndex
1946 for (uint8_t nIndex = 0; nIndex < channelParams->m_reducedClusterNumber; nIndex++)
1947 {
1948 for (uint8_t mIndex = 0; mIndex < table3gpp->m_raysPerCluster; mIndex++)
1949 {
1950 DoubleVector initialPhase = channelParams->m_clusterPhase[nIndex][mIndex];
1951 NS_ASSERT(4 <= initialPhase.size());
1952 double k = channelParams->m_crossPolarizationPowerRatios[nIndex][mIndex];
1953
1954 // cache the component of the "rays" terms which depend on the random angle of arrivals
1955 // and departures and initial phases only
1956 auto [rxFieldPatternPhi, rxFieldPatternTheta] = uAntenna->GetElementFieldPattern(
1957 Angles(channelParams->m_rayAoaRadian[nIndex][mIndex],
1958 channelParams->m_rayZoaRadian[nIndex][mIndex]));
1959 auto [txFieldPatternPhi, txFieldPatternTheta] = sAntenna->GetElementFieldPattern(
1960 Angles(channelParams->m_rayAodRadian[nIndex][mIndex],
1961 channelParams->m_rayZodRadian[nIndex][mIndex]));
1962 raysPreComp(nIndex, mIndex) =
1963 std::complex<double>(cos(initialPhase[0]), sin(initialPhase[0])) *
1964 rxFieldPatternTheta * txFieldPatternTheta +
1965 std::complex<double>(cos(initialPhase[1]), sin(initialPhase[1])) *
1966 std::sqrt(1.0 / k) * rxFieldPatternTheta * txFieldPatternPhi +
1967 std::complex<double>(cos(initialPhase[2]), sin(initialPhase[2])) *
1968 std::sqrt(1.0 / k) * rxFieldPatternPhi * txFieldPatternTheta +
1969 std::complex<double>(cos(initialPhase[3]), sin(initialPhase[3])) *
1970 rxFieldPatternPhi * txFieldPatternPhi;
1971
1972 // cache the component of the "rxPhaseDiff" terms which depend on the random angle of
1973 // arrivals only
1974 double sinRayZoa = sin(rayZoaRadian[nIndex][mIndex]);
1975 double sinRayAoa = cos(rayAoaRadian[nIndex][mIndex]);
1976 double cosRayAoa = cos(rayAoaRadian[nIndex][mIndex]);
1977 sinCosA[nIndex][mIndex] = sinRayZoa * cosRayAoa;
1978 sinSinA[nIndex][mIndex] = sinRayZoa * sinRayAoa;
1979 cosZoA[nIndex][mIndex] = cos(rayZoaRadian[nIndex][mIndex]);
1980
1981 // cache the component of the "txPhaseDiff" terms which depend on the random angle of
1982 // departure only
1983 double sinRayZod = sin(rayZodRadian[nIndex][mIndex]);
1984 double sinRayAod = cos(rayAodRadian[nIndex][mIndex]);
1985 double cosRayAod = cos(rayAodRadian[nIndex][mIndex]);
1986 sinCosD[nIndex][mIndex] = sinRayZod * cosRayAod;
1987 sinSinD[nIndex][mIndex] = sinRayZod * sinRayAod;
1988 cosZoD[nIndex][mIndex] = cos(rayZodRadian[nIndex][mIndex]);
1989 }
1990 }
1991
1992 // The following for loops computes the channel coefficients
1993 // Keeps track of how many sub-clusters have been added up to now
1994 uint8_t numSubClustersAdded = 0;
1995 for (uint8_t nIndex = 0; nIndex < channelParams->m_reducedClusterNumber; nIndex++)
1996 {
1997 for (size_t uIndex = 0; uIndex < uSize; uIndex++)
1998 {
1999 Vector uLoc = uAntenna->GetElementLocation(uIndex);
2000
2001 for (size_t sIndex = 0; sIndex < sSize; sIndex++)
2002 {
2003 Vector sLoc = sAntenna->GetElementLocation(sIndex);
2004 // Compute the N-2 weakest cluster, assuming 0 slant angle and a
2005 // polarization slant angle configured in the array (7.5-22)
2006 if (nIndex != channelParams->m_cluster1st && nIndex != channelParams->m_cluster2nd)
2007 {
2008 std::complex<double> rays(0, 0);
2009 for (uint8_t mIndex = 0; mIndex < table3gpp->m_raysPerCluster; mIndex++)
2010 {
2011 // lambda_0 is accounted in the antenna spacing uLoc and sLoc.
2012 double rxPhaseDiff =
2013 2 * M_PI *
2014 (sinCosA[nIndex][mIndex] * uLoc.x + sinSinA[nIndex][mIndex] * uLoc.y +
2015 cosZoA[nIndex][mIndex] * uLoc.z);
2016
2017 double txPhaseDiff =
2018 2 * M_PI *
2019 (sinCosD[nIndex][mIndex] * sLoc.x + sinSinD[nIndex][mIndex] * sLoc.y +
2020 cosZoD[nIndex][mIndex] * sLoc.z);
2021
2022 // NOTE Doppler is computed in the CalcBeamformingGain function and is
2023 // simplified to only account for the center angle of each cluster.
2024 rays += raysPreComp(nIndex, mIndex) *
2025 std::complex<double>(cos(rxPhaseDiff), sin(rxPhaseDiff)) *
2026 std::complex<double>(cos(txPhaseDiff), sin(txPhaseDiff));
2027 }
2028 rays *=
2029 sqrt(channelParams->m_clusterPower[nIndex] / table3gpp->m_raysPerCluster);
2030 hUsn(uIndex, sIndex, nIndex) = rays;
2031 }
2032 else //(7.5-28)
2033 {
2034 std::complex<double> raysSub1(0, 0);
2035 std::complex<double> raysSub2(0, 0);
2036 std::complex<double> raysSub3(0, 0);
2037
2038 for (uint8_t mIndex = 0; mIndex < table3gpp->m_raysPerCluster; mIndex++)
2039 {
2040 // ZML:Just remind me that the angle offsets for the 3 subclusters were not
2041 // generated correctly.
2042 double rxPhaseDiff =
2043 2 * M_PI *
2044 (sinCosA[nIndex][mIndex] * uLoc.x + sinSinA[nIndex][mIndex] * uLoc.y +
2045 cosZoA[nIndex][mIndex] * uLoc.z);
2046
2047 double txPhaseDiff =
2048 2 * M_PI *
2049 (sinCosD[nIndex][mIndex] * sLoc.x + sinSinD[nIndex][mIndex] * sLoc.y +
2050 cosZoD[nIndex][mIndex] * sLoc.z);
2051
2052 std::complex<double> raySub =
2053 raysPreComp(nIndex, mIndex) *
2054 std::complex<double>(cos(rxPhaseDiff), sin(rxPhaseDiff)) *
2055 std::complex<double>(cos(txPhaseDiff), sin(txPhaseDiff));
2056
2057 switch (mIndex)
2058 {
2059 case 9:
2060 case 10:
2061 case 11:
2062 case 12:
2063 case 17:
2064 case 18:
2065 raysSub2 += raySub;
2066 break;
2067 case 13:
2068 case 14:
2069 case 15:
2070 case 16:
2071 raysSub3 += raySub;
2072 break;
2073 default: // case 1,2,3,4,5,6,7,8,19,20
2074 raysSub1 += raySub;
2075 break;
2076 }
2077 }
2078 raysSub1 *=
2079 sqrt(channelParams->m_clusterPower[nIndex] / table3gpp->m_raysPerCluster);
2080 raysSub2 *=
2081 sqrt(channelParams->m_clusterPower[nIndex] / table3gpp->m_raysPerCluster);
2082 raysSub3 *=
2083 sqrt(channelParams->m_clusterPower[nIndex] / table3gpp->m_raysPerCluster);
2084 hUsn(uIndex, sIndex, nIndex) = raysSub1;
2085 hUsn(uIndex,
2086 sIndex,
2087 channelParams->m_reducedClusterNumber + numSubClustersAdded) = raysSub2;
2088 hUsn(uIndex,
2089 sIndex,
2090 channelParams->m_reducedClusterNumber + numSubClustersAdded + 1) =
2091 raysSub3;
2092 }
2093 }
2094 }
2095 if (nIndex == channelParams->m_cluster1st || nIndex == channelParams->m_cluster2nd)
2096 {
2097 numSubClustersAdded += 2;
2098 }
2099 }
2100
2101 if (channelParams->m_losCondition == ChannelCondition::LOS) //(7.5-29) && (7.5-30)
2102 {
2103 double lambda = 3.0e8 / m_frequency; // the wavelength of the carrier frequency
2104 std::complex<double> phaseDiffDueToDistance(cos(-2 * M_PI * distance3D / lambda),
2105 sin(-2 * M_PI * distance3D / lambda));
2106
2107 const double sinUAngleIncl = sin(uAngle.GetInclination());
2108 const double cosUAngleIncl = cos(uAngle.GetInclination());
2109 const double sinUAngleAz = sin(uAngle.GetAzimuth());
2110 const double cosUAngleAz = cos(uAngle.GetAzimuth());
2111 const double sinSAngleIncl = sin(sAngle.GetInclination());
2112 const double cosSAngleIncl = cos(sAngle.GetInclination());
2113 const double sinSAngleAz = sin(sAngle.GetAzimuth());
2114 const double cosSAngleAz = cos(sAngle.GetAzimuth());
2115
2116 for (size_t uIndex = 0; uIndex < uSize; uIndex++)
2117 {
2118 Vector uLoc = uAntenna->GetElementLocation(uIndex);
2119 double rxPhaseDiff = 2 * M_PI *
2120 (sinUAngleIncl * cosUAngleAz * uLoc.x +
2121 sinUAngleIncl * sinUAngleAz * uLoc.y + cosUAngleIncl * uLoc.z);
2122
2123 for (size_t sIndex = 0; sIndex < sSize; sIndex++)
2124 {
2125 Vector sLoc = sAntenna->GetElementLocation(sIndex);
2126 std::complex<double> ray(0, 0);
2127 double txPhaseDiff =
2128 2 * M_PI *
2129 (sinSAngleIncl * cosSAngleAz * sLoc.x + sinSAngleIncl * sinSAngleAz * sLoc.y +
2130 cosSAngleIncl * sLoc.z);
2131
2132 auto [rxFieldPatternPhi, rxFieldPatternTheta] = uAntenna->GetElementFieldPattern(
2133 Angles(uAngle.GetAzimuth(), uAngle.GetInclination()));
2134 auto [txFieldPatternPhi, txFieldPatternTheta] = sAntenna->GetElementFieldPattern(
2135 Angles(sAngle.GetAzimuth(), sAngle.GetInclination()));
2136
2137 ray = (rxFieldPatternTheta * txFieldPatternTheta -
2138 rxFieldPatternPhi * txFieldPatternPhi) *
2139 phaseDiffDueToDistance *
2140 std::complex<double>(cos(rxPhaseDiff), sin(rxPhaseDiff)) *
2141 std::complex<double>(cos(txPhaseDiff), sin(txPhaseDiff));
2142
2143 double kLinear = pow(10, channelParams->m_K_factor / 10.0);
2144 // the LOS path should be attenuated if blockage is enabled.
2145 hUsn(uIndex, sIndex, 0) =
2146 sqrt(1.0 / (kLinear + 1)) * hUsn(uIndex, sIndex, 0) +
2147 sqrt(kLinear / (1 + kLinear)) * ray /
2148 pow(10,
2149 channelParams->m_attenuation_dB[0] / 10.0); //(7.5-30) for tau = tau1
2150 for (uint16_t nIndex = 1; nIndex < hUsn.GetNumPages(); nIndex++)
2151 {
2152 hUsn(uIndex, sIndex, nIndex) *=
2153 sqrt(1.0 / (kLinear + 1)); //(7.5-30) for tau = tau2...tauN
2154 }
2155 }
2156 }
2157 }
2158
2159 NS_LOG_DEBUG("Husn (sAntenna, uAntenna):" << sAntenna->GetId() << ", " << uAntenna->GetId());
2160 for (uint16_t cIndex = 0; cIndex < hUsn.GetNumPages(); cIndex++)
2161 {
2162 for (uint16_t rowIdx = 0; rowIdx < hUsn.GetNumRows(); rowIdx++)
2163 {
2164 for (uint16_t colIdx = 0; colIdx < hUsn.GetNumCols(); colIdx++)
2165 {
2166 NS_LOG_DEBUG(" " << hUsn(rowIdx, colIdx, cIndex) << ",");
2167 }
2168 }
2169 }
2170
2171 NS_LOG_INFO("size of coefficient matrix (rows, columns, clusters) = ("
2172 << hUsn.GetNumRows() << ", " << hUsn.GetNumCols() << ", " << hUsn.GetNumPages()
2173 << ")");
2174 channelMatrix->m_channel = hUsn;
2175 return channelMatrix;
2176}
2177
2178std::pair<double, double>
2179ThreeGppChannelModel::WrapAngles(double azimuthRad, double inclinationRad)
2180{
2181 inclinationRad = WrapTo2Pi(inclinationRad);
2182 if (inclinationRad > M_PI)
2183 {
2184 // inclination must be in [0, M_PI]
2185 inclinationRad -= M_PI;
2186 azimuthRad += M_PI;
2187 }
2188
2189 azimuthRad = WrapTo2Pi(azimuthRad);
2190
2191 NS_ASSERT_MSG(0 <= inclinationRad && inclinationRad <= M_PI,
2192 "inclinationRad=" << inclinationRad << " not valid, should be in [0, pi]");
2193 NS_ASSERT_MSG(0 <= azimuthRad && azimuthRad <= 2 * M_PI,
2194 "azimuthRad=" << azimuthRad << " not valid, should be in [0, 2*pi]");
2195
2196 return std::make_pair(azimuthRad, inclinationRad);
2197}
2198
2202 const DoubleVector& clusterAOA,
2203 const DoubleVector& clusterZOA) const
2204{
2205 NS_LOG_FUNCTION(this);
2206
2207 DoubleVector powerAttenuation;
2208 uint8_t clusterNum = clusterAOA.size();
2209 for (uint8_t cInd = 0; cInd < clusterNum; cInd++)
2210 {
2211 powerAttenuation.push_back(0); // Initial power attenuation for all clusters to be 0 dB;
2212 }
2213 // step a: the number of non-self blocking blockers is stored in m_numNonSelfBlocking.
2214
2215 // step b:Generate the size and location of each blocker
2216 // generate self blocking (i.e., for blockage from the human body)
2217 // table 7.6.4.1-1 Self-blocking region parameters.
2218 // Defaults: landscape mode
2219 double phiSb = 40;
2220 double xSb = 160;
2221 double thetaSb = 110;
2222 double ySb = 75;
2223 if (m_portraitMode)
2224 {
2225 phiSb = 260;
2226 xSb = 120;
2227 thetaSb = 100;
2228 ySb = 80;
2229 }
2230
2231 // generate or update non-self blocking
2232 if (channelParams->m_nonSelfBlocking.empty()) // generate new blocking regions
2233 {
2234 for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
2235 {
2236 // draw value from table 7.6.4.1-2 Blocking region parameters
2237 DoubleVector table;
2238 table.push_back(m_normalRv->GetValue()); // phi_k: store the normal RV that will be
2239 // mapped to uniform (0,360) later.
2240 if (m_scenario == "InH-OfficeMixed" || m_scenario == "InH-OfficeOpen")
2241 {
2242 table.push_back(m_uniformRv->GetValue(15, 45)); // x_k
2243 table.push_back(90); // Theta_k
2244 table.push_back(m_uniformRv->GetValue(5, 15)); // y_k
2245 table.push_back(2); // r
2246 }
2247 else
2248 {
2249 table.push_back(m_uniformRv->GetValue(5, 15)); // x_k
2250 table.push_back(90); // Theta_k
2251 table.push_back(5); // y_k
2252 table.push_back(10); // r
2253 }
2254 channelParams->m_nonSelfBlocking.push_back(table);
2255 }
2256 }
2257 else
2258 {
2259 double deltaX = sqrt(pow(channelParams->m_preLocUT.x - channelParams->m_locUT.x, 2) +
2260 pow(channelParams->m_preLocUT.y - channelParams->m_locUT.y, 2));
2261 // if deltaX and speed are both 0, the autocorrelation is 1, skip updating
2262 if (deltaX > 1e-6 || m_blockerSpeed > 1e-6)
2263 {
2264 double corrDis;
2265 // draw value from table 7.6.4.1-4: Spatial correlation distance for different
2266 // m_scenarios.
2267 if (m_scenario == "InH-OfficeMixed" || m_scenario == "InH-OfficeOpen")
2268 {
2269 // InH, correlation distance = 5;
2270 corrDis = 5;
2271 }
2272 else
2273 {
2274 if (channelParams->m_o2iCondition == ChannelCondition::O2I) // outdoor to indoor
2275 {
2276 corrDis = 5;
2277 }
2278 else // LOS or NLOS
2279 {
2280 corrDis = 10;
2281 }
2282 }
2283 double R;
2284 if (m_blockerSpeed > 1e-6) // speed not equal to 0
2285 {
2286 double corrT = corrDis / m_blockerSpeed;
2287 R = exp(-1 * (deltaX / corrDis +
2288 (Now().GetSeconds() - channelParams->m_generatedTime.GetSeconds()) /
2289 corrT));
2290 }
2291 else
2292 {
2293 R = exp(-1 * (deltaX / corrDis));
2294 }
2295
2296 NS_LOG_INFO("Distance change:"
2297 << deltaX << " Speed:" << m_blockerSpeed << " Time difference:"
2298 << Now().GetSeconds() - channelParams->m_generatedTime.GetSeconds()
2299 << " correlation:" << R);
2300
2301 // In order to generate correlated uniform random variables, we first generate
2302 // correlated normal random variables and map the normal RV to uniform RV. Notice the
2303 // correlation will change if the RV is transformed from normal to uniform. To
2304 // compensate the distortion, the correlation of the normal RV is computed such that the
2305 // uniform RV would have the desired correlation when transformed from normal RV.
2306
2307 // The following formula was obtained from MATLAB numerical simulation.
2308
2309 if (R * R * (-0.069) + R * 1.074 - 0.002 <
2310 1) // transform only when the correlation of normal RV is smaller than 1
2311 {
2312 R = R * R * (-0.069) + R * 1.074 - 0.002;
2313 }
2314 for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
2315 {
2316 // Generate a new correlated normal RV with the following formula
2317 channelParams->m_nonSelfBlocking[blockInd][PHI_INDEX] =
2318 R * channelParams->m_nonSelfBlocking[blockInd][PHI_INDEX] +
2319 sqrt(1 - R * R) * m_normalRv->GetValue();
2320 }
2321 }
2322 }
2323
2324 // step c: Determine the attenuation of each blocker due to blockers
2325 for (uint8_t cInd = 0; cInd < clusterNum; cInd++)
2326 {
2327 NS_ASSERT_MSG(clusterAOA[cInd] >= 0 && clusterAOA[cInd] <= 360,
2328 "the AOA should be the range of [0,360]");
2329 NS_ASSERT_MSG(clusterZOA[cInd] >= 0 && clusterZOA[cInd] <= 180,
2330 "the ZOA should be the range of [0,180]");
2331
2332 // check self blocking
2333 NS_LOG_INFO("AOA=" << clusterAOA[cInd] << " Block Region[" << phiSb - xSb / 2.0 << ","
2334 << phiSb + xSb / 2.0 << "]");
2335 NS_LOG_INFO("ZOA=" << clusterZOA[cInd] << " Block Region[" << thetaSb - ySb / 2.0 << ","
2336 << thetaSb + ySb / 2.0 << "]");
2337 if (std::abs(clusterAOA[cInd] - phiSb) < (xSb / 2.0) &&
2338 std::abs(clusterZOA[cInd] - thetaSb) < (ySb / 2.0))
2339 {
2340 powerAttenuation[cInd] += 30; // attenuate by 30 dB.
2341 NS_LOG_INFO("Cluster[" << +cInd
2342 << "] is blocked by self blocking region and reduce 30 dB power,"
2343 "the attenuation is ["
2344 << powerAttenuation[cInd] << " dB]");
2345 }
2346
2347 // check non-self blocking
2348 for (uint16_t blockInd = 0; blockInd < m_numNonSelfBlocking; blockInd++)
2349 {
2350 // The normal RV is transformed to uniform RV with the desired correlation.
2351 double phiK =
2352 (0.5 * erfc(-1 * channelParams->m_nonSelfBlocking[blockInd][PHI_INDEX] / sqrt(2))) *
2353 360;
2354 while (phiK > 360)
2355 {
2356 phiK -= 360;
2357 }
2358
2359 while (phiK < 0)
2360 {
2361 phiK += 360;
2362 }
2363
2364 double xK = channelParams->m_nonSelfBlocking[blockInd][X_INDEX];
2365 double thetaK = channelParams->m_nonSelfBlocking[blockInd][THETA_INDEX];
2366 double yK = channelParams->m_nonSelfBlocking[blockInd][Y_INDEX];
2367
2368 NS_LOG_INFO("AOA=" << clusterAOA[cInd] << " Block Region[" << phiK - xK << ","
2369 << phiK + xK << "]");
2370 NS_LOG_INFO("ZOA=" << clusterZOA[cInd] << " Block Region[" << thetaK - yK << ","
2371 << thetaK + yK << "]");
2372
2373 if (std::abs(clusterAOA[cInd] - phiK) < (xK) &&
2374 std::abs(clusterZOA[cInd] - thetaK) < (yK))
2375 {
2376 double A1 = clusterAOA[cInd] - (phiK + xK / 2.0); //(7.6-24)
2377 double A2 = clusterAOA[cInd] - (phiK - xK / 2.0); //(7.6-25)
2378 double Z1 = clusterZOA[cInd] - (thetaK + yK / 2.0); //(7.6-26)
2379 double Z2 = clusterZOA[cInd] - (thetaK - yK / 2.0); //(7.6-27)
2380 int signA1;
2381 int signA2;
2382 int signZ1;
2383 int signZ2;
2384 // draw sign for the above parameters according to table 7.6.4.1-3 Description of
2385 // signs
2386 if (xK / 2.0 < clusterAOA[cInd] - phiK && clusterAOA[cInd] - phiK <= xK)
2387 {
2388 signA1 = -1;
2389 }
2390 else
2391 {
2392 signA1 = 1;
2393 }
2394 if (-1 * xK < clusterAOA[cInd] - phiK && clusterAOA[cInd] - phiK <= -1 * xK / 2.0)
2395 {
2396 signA2 = -1;
2397 }
2398 else
2399 {
2400 signA2 = 1;
2401 }
2402
2403 if (yK / 2.0 < clusterZOA[cInd] - thetaK && clusterZOA[cInd] - thetaK <= yK)
2404 {
2405 signZ1 = -1;
2406 }
2407 else
2408 {
2409 signZ1 = 1;
2410 }
2411 if (-1 * yK < clusterZOA[cInd] - thetaK &&
2412 clusterZOA[cInd] - thetaK <= -1 * yK / 2.0)
2413 {
2414 signZ2 = -1;
2415 }
2416 else
2417 {
2418 signZ2 = 1;
2419 }
2420 double lambda = 3e8 / m_frequency;
2421 double fA1 =
2422 atan(signA1 * M_PI / 2.0 *
2423 sqrt(M_PI / lambda * channelParams->m_nonSelfBlocking[blockInd][R_INDEX] *
2424 (1.0 / cos(DegreesToRadians(A1)) - 1))) /
2425 M_PI; //(7.6-23)
2426 double fA2 =
2427 atan(signA2 * M_PI / 2.0 *
2428 sqrt(M_PI / lambda * channelParams->m_nonSelfBlocking[blockInd][R_INDEX] *
2429 (1.0 / cos(DegreesToRadians(A2)) - 1))) /
2430 M_PI;
2431 double fZ1 =
2432 atan(signZ1 * M_PI / 2.0 *
2433 sqrt(M_PI / lambda * channelParams->m_nonSelfBlocking[blockInd][R_INDEX] *
2434 (1.0 / cos(DegreesToRadians(Z1)) - 1))) /
2435 M_PI;
2436 double fZ2 =
2437 atan(signZ2 * M_PI / 2.0 *
2438 sqrt(M_PI / lambda * channelParams->m_nonSelfBlocking[blockInd][R_INDEX] *
2439 (1.0 / cos(DegreesToRadians(Z2)) - 1))) /
2440 M_PI;
2441 double lDb = -20 * log10(1 - (fA1 + fA2) * (fZ1 + fZ2)); //(7.6-22)
2442 powerAttenuation[cInd] += lDb;
2443 NS_LOG_INFO("Cluster[" << +cInd << "] is blocked by no-self blocking, the loss is ["
2444 << lDb << "] dB");
2445 }
2446 }
2447 }
2448 return powerAttenuation;
2449}
2450
2451void
2452ThreeGppChannelModel::Shuffle(double* first, double* last) const
2453{
2454 for (auto i = (last - first) - 1; i > 0; --i)
2455 {
2456 std::swap(first[i], first[m_uniformRvShuffle->GetInteger(0, i)]);
2457 }
2458}
2459
2460int64_t
2462{
2463 NS_LOG_FUNCTION(this << stream);
2464 m_normalRv->SetStream(stream);
2465 m_uniformRv->SetStream(stream + 1);
2466 m_uniformRvShuffle->SetStream(stream + 2);
2467 m_uniformRvDoppler->SetStream(stream + 3);
2468 return 4;
2469}
2470
2471} // namespace ns3
#define min(a, b)
Definition: 80211b.c:42
double f(double x, void *params)
Definition: 80211b.c:71
#define max(a, b)
Definition: 80211b.c:43
Class holding the azimuth and inclination angles of spherical coordinates.
Definition: angles.h:118
double GetInclination() const
Getter for inclination angle.
Definition: angles.cc:216
double GetAzimuth() const
Getter for azimuth angle.
Definition: angles.cc:210
AttributeValue implementation for Boolean.
Definition: boolean.h:37
This class can be used to hold variables of floating point type such as 'double' or 'float'.
Definition: double.h:42
Hold a signed integer type.
Definition: integer.h:45
MatrixArray class inherits ValArray class and provides additional interfaces to ValArray which enable...
Definition: matrix-array.h:83
This is an interface for a channel model that can be described by a channel matrix,...
std::vector< double > DoubleVector
Type definition for vectors of doubles.
std::vector< Double2DVector > Double3DVector
Type definition for 3D matrices of doubles.
std::vector< DoubleVector > Double2DVector
Type definition for matrices of doubles.
static uint64_t GetKey(uint32_t a, uint32_t b)
Generate a unique value for the pair of unsigned integer of 32 bits, where the order does not matter,...
A network Node.
Definition: node.h:56
uint32_t GetId() const
Definition: node.cc:117
Hold objects of type Ptr<T>.
Definition: pointer.h:37
Smart pointer class similar to boost::intrusive_ptr.
Definition: ptr.h:78
void SetStream(int64_t stream)
Specifies the stream number for the RngStream.
static Time Now()
Return the current simulation virtual time.
Definition: simulator.cc:199
Hold variables of type string.
Definition: string.h:56
DoubleVector CalcAttenuationOfBlockage(const Ptr< ThreeGppChannelModel::ThreeGppChannelParams > channelParams, const DoubleVector &clusterAOA, const DoubleVector &clusterZOA) const
Applies the blockage model A described in 3GPP TR 38.901.
int64_t AssignStreams(int64_t stream)
Assign a fixed random variable stream number to the random variables used by this model.
bool m_portraitMode
true if portrait mode, false if landscape
bool ChannelParamsNeedsUpdate(Ptr< const ThreeGppChannelParams > channelParams, Ptr< const ChannelCondition > channelCondition) const
Check if the channel params has to be updated.
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.
Ptr< NormalRandomVariable > m_normalRv
normal random variable
static const uint8_t Y_INDEX
index of the Y value in the m_nonSelfBlocking array
bool m_blockage
enables the blockage model A
Ptr< const ChannelParams > GetParams(Ptr< const MobilityModel > aMob, Ptr< const MobilityModel > bMob) const override
Looks for the channel params associated to the aMob and bMob pair in m_channelParamsMap.
~ThreeGppChannelModel() override
Destructor.
bool ChannelMatrixNeedsUpdate(Ptr< const ThreeGppChannelParams > channelParams, Ptr< const ChannelMatrix > channelMatrix)
Check if the channel matrix has to be updated (it needs update when the channel params generation tim...
static const uint8_t THETA_INDEX
index of the THETA value in the m_nonSelfBlocking array
std::unordered_map< uint64_t, Ptr< ThreeGppChannelParams > > m_channelParamsMap
map containing the common channel parameters per pair of nodes, the key of this map is reciprocal and...
static std::pair< double, double > WrapAngles(double azimuthRad, double inclinationRad)
Wrap an (azimuth, inclination) angle pair in a valid range.
double m_blockerSpeed
the blocker speed
Ptr< const ChannelMatrix > GetChannel(Ptr< const MobilityModel > aMob, Ptr< const MobilityModel > bMob, Ptr< const PhasedArrayModel > aAntenna, Ptr< const PhasedArrayModel > bAntenna) override
Looks for the channel matrix associated to the aMob and bMob pair in m_channelMatrixMap.
void SetFrequency(double f)
Sets the center frequency of the model.
std::unordered_map< uint64_t, Ptr< ChannelMatrix > > m_channelMatrixMap
map containing the channel realizations per pair of PhasedAntennaArray instances, the key of this map...
Ptr< UniformRandomVariable > m_uniformRv
uniform random variable
void DoDispose() override
Destructor implementation.
void SetScenario(const std::string &scenario)
Sets the propagation scenario.
void SetChannelConditionModel(Ptr< ChannelConditionModel > model)
Set the channel condition model.
Ptr< UniformRandomVariable > m_uniformRvDoppler
uniform random variable, used to compute the additional Doppler contribution
uint16_t m_numNonSelfBlocking
number of non-self-blocking regions
std::string GetScenario() const
Returns the propagation scenario.
virtual Ptr< ChannelMatrix > GetNewChannel(Ptr< const ThreeGppChannelParams > channelParams, Ptr< const ParamsTable > table3gpp, const Ptr< const MobilityModel > sMob, const Ptr< const MobilityModel > uMob, Ptr< const PhasedArrayModel > sAntenna, Ptr< const PhasedArrayModel > uAntenna) const
Compute the channel matrix between two nodes a and b, and their antenna arrays aAntenna and bAntenna ...
static const uint8_t PHI_INDEX
index of the PHI value in the m_nonSelfBlocking array
double m_frequency
the operating frequency
double m_vScatt
value used to compute the additional Doppler contribution for the delayed paths
Ptr< ChannelConditionModel > GetChannelConditionModel() const
Get the associated channel condition model.
Ptr< ChannelConditionModel > m_channelConditionModel
the channel condition model
std::string m_scenario
the 3GPP scenario
static const uint8_t R_INDEX
index of the R value in the m_nonSelfBlocking array
static TypeId GetTypeId()
Get the type ID.
void Shuffle(double *first, double *last) const
Shuffle the elements of a simple sequence container of type double.
Ptr< ThreeGppChannelParams > GenerateChannelParameters(const Ptr< const ChannelCondition > channelCondition, const Ptr< const ParamsTable > table3gpp, const Ptr< const MobilityModel > aMob, const Ptr< const MobilityModel > bMob) const
Prepare 3gpp channel parameters among the nodes a and b.
double GetFrequency() const
Returns the center frequency.
Time m_updatePeriod
the channel update period
static const uint8_t X_INDEX
index of the X value in the m_nonSelfBlocking array
Ptr< UniformRandomVariable > m_uniformRvShuffle
uniform random variable used to shuffle array in GetNewChannel
@ NS
nanosecond
Definition: nstime.h:119
bool IsZero() const
Exactly equivalent to t == 0.
Definition: nstime.h:314
AttributeValue implementation for Time.
Definition: nstime.h:1423
a unique identifier for an interface.
Definition: type-id.h:60
TypeId SetGroupName(std::string groupName)
Set the group name.
Definition: type-id.cc:943
TypeId SetParent(TypeId tid)
Set the parent TypeId.
Definition: type-id.cc:935
double GetValue(double min, double max)
Get the next random value drawn from the distribution.
uint32_t GetInteger(uint32_t min, uint32_t max)
Get the next random value drawn from the distribution.
#define NS_ASSERT(condition)
At runtime, in debugging builds, if this condition is not true, the program prints the source file,...
Definition: assert.h:66
#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:86
Ptr< const AttributeAccessor > MakeBooleanAccessor(T1 a1)
Definition: boolean.h:86
Ptr< const AttributeChecker > MakeBooleanChecker()
Definition: boolean.cc:124
Ptr< const AttributeAccessor > MakeDoubleAccessor(T1 a1)
Definition: double.h:43
Ptr< const AttributeAccessor > MakeIntegerAccessor(T1 a1)
Definition: integer.h:46
Ptr< const AttributeAccessor > MakePointerAccessor(T1 a1)
Definition: pointer.h:231
Ptr< const AttributeChecker > MakeStringChecker()
Definition: string.cc:30
Ptr< const AttributeAccessor > MakeStringAccessor(T1 a1)
Definition: string.h:57
Ptr< const AttributeChecker > MakeTimeChecker()
Helper to make an unbounded Time checker.
Definition: nstime.h:1444
Ptr< const AttributeAccessor > MakeTimeAccessor(T1 a1)
Definition: nstime.h:1424
#define NS_FATAL_ERROR(msg)
Report a fatal error with a message and terminate.
Definition: fatal-error.h:179
#define NS_LOG_COMPONENT_DEFINE(name)
Define a Log component with a specific name.
Definition: log.h:202
#define NS_LOG_DEBUG(msg)
Use NS_LOG to output a message of level LOG_DEBUG.
Definition: log.h:268
#define NS_LOG_FUNCTION(parameters)
If log level LOG_FUNCTION is enabled, this macro will output all input parameters separated by ",...
#define NS_LOG_WARN(msg)
Use NS_LOG to output a message of level LOG_WARN.
Definition: log.h:261
#define NS_LOG_INFO(msg)
Use NS_LOG to output a message of level LOG_INFO.
Definition: log.h:275
#define NS_OBJECT_ENSURE_REGISTERED(type)
Register an Object subclass with the TypeId system.
Definition: object-base.h:46
Time Now()
create an ns3::Time instance which contains the current simulation time.
Definition: simulator.cc:296
Time MilliSeconds(uint64_t value)
Construct a Time in the indicated unit.
Definition: nstime.h:1348
Definition: first.py:1
Every class exported by the ns3 library is enclosed in the ns3 namespace.
static const double offSetAlpha[20]
The ray offset angles within a cluster, given for rms angle spread normalized to 1.
static const double sqrtC_RMa_O2I[6][6]
The square root matrix for RMa O2I, which is generated using the Cholesky decomposition according to ...
static const double sqrtC_UMi_LOS[7][7]
The square root matrix for UMi LOS, which is generated using the Cholesky decomposition according to ...
static const double sqrtC_office_LOS[7][7]
The square root matrix for Indoor-Office LOS, which is generated using the Cholesky decomposition acc...
static const double sqrtC_UMa_O2I[6][6]
The square root matrix for UMa O2I, which is generated using the Cholesky decomposition according to ...
static const double sqrtC_RMa_NLOS[6][6]
The square root matrix for RMa NLOS, which is generated using the Cholesky decomposition according to...
static const double sqrtC_UMa_LOS[7][7]
The square root matrix for UMa LOS, which is generated using the Cholesky decomposition according to ...
static const double sqrtC_UMi_NLOS[6][6]
The square root matrix for UMi NLOS, which is generated using the Cholesky decomposition according to...
static const double sqrtC_RMa_LOS[7][7]
The square root matrix for RMa LOS, which is generated using the Cholesky decomposition according to ...
double DegreesToRadians(double degrees)
converts degrees to radians
Definition: angles.cc:39
static const double sqrtC_UMi_O2I[6][6]
The square root matrix for UMi O2I, which is generated using the Cholesky decomposition according to ...
static const double sqrtC_office_NLOS[6][6]
The square root matrix for Indoor-Office NLOS, which is generated using the Cholesky decomposition ac...
static const double sqrtC_UMa_NLOS[6][6]
The square root matrix for UMa NLOS, which is generated using the Cholesky decomposition according to...
double WrapTo2Pi(double a)
Wrap angle in [0, 2*M_PI)
Definition: angles.cc:102
double RadiansToDegrees(double radians)
converts radians to degrees
Definition: angles.cc:45