379 lines
11 KiB
Go
379 lines
11 KiB
Go
|
package simulation
|
|||
|
|
|||
|
import (
|
|||
|
"encoding/json"
|
|||
|
"fmt"
|
|||
|
"math"
|
|||
|
"math/rand"
|
|||
|
"sync"
|
|||
|
|
|||
|
vector "moxitech/drone/simulation-go/pkg/vector"
|
|||
|
"moxitech/drone/simulation-go/utils"
|
|||
|
|
|||
|
"github.com/dhconnelly/rtreego" // R-Tree для пространственного индекса
|
|||
|
)
|
|||
|
|
|||
|
type UAV struct {
|
|||
|
Path map[int][]float64
|
|||
|
Speed *vector.Vector3
|
|||
|
Freq float64 // Частота передачи
|
|||
|
Radius float64 // Радиус покрытия в километрах
|
|||
|
}
|
|||
|
|
|||
|
type CGS struct {
|
|||
|
Position *vector.Vector3
|
|||
|
Radius float64 // Радиус покрытия в километрах
|
|||
|
Frequency float64 // Частота передачи в МГц
|
|||
|
}
|
|||
|
|
|||
|
type Simulation struct {
|
|||
|
elevationDefault float64
|
|||
|
distanceDefault float64
|
|||
|
uavs map[int]*UAV
|
|||
|
cgs map[int]*CGS
|
|||
|
heights map[string]float64 // высоты в формате "lat,lon": высота
|
|||
|
idGen func() int
|
|||
|
timeline []int
|
|||
|
mutex sync.Mutex
|
|||
|
rtree *rtreego.Rtree
|
|||
|
results []SimulationResult
|
|||
|
}
|
|||
|
|
|||
|
type SimulationResult struct {
|
|||
|
Time int `json:"time"`
|
|||
|
UAVs map[int]UAVResult `json:"uavs"`
|
|||
|
CGSs map[int]CGSResult `json:"cgs"`
|
|||
|
}
|
|||
|
|
|||
|
type UAVResult struct {
|
|||
|
Position *vector.Vector3 `json:"position"`
|
|||
|
Freq float64 `json:"frequency"`
|
|||
|
Radius float64 `json:"radius"`
|
|||
|
DataRate float64 `json:"data_rate"`
|
|||
|
}
|
|||
|
|
|||
|
type CGSResult struct {
|
|||
|
Position *vector.Vector3 `json:"position"`
|
|||
|
Radius float64 `json:"radius"`
|
|||
|
Frequency float64 `json:"frequency"`
|
|||
|
}
|
|||
|
|
|||
|
func NewSimulation() *Simulation {
|
|||
|
return &Simulation{
|
|||
|
elevationDefault: 85.0 * (math.Pi / 180), // 15 degrees in radians
|
|||
|
distanceDefault: 1000.0, // km
|
|||
|
uavs: make(map[int]*UAV),
|
|||
|
cgs: make(map[int]*CGS),
|
|||
|
heights: make(map[string]float64),
|
|||
|
idGen: utils.GetIDGenerator(),
|
|||
|
rtree: rtreego.NewTree(2, 25, 50), // Инициализация R-Tree
|
|||
|
results: []SimulationResult{},
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
// Adding CGS (ground station) with coverage radius and frequency
|
|||
|
func (s *Simulation) AddCGS(pos []float64, radius float64, frequency float64) {
|
|||
|
s.mutex.Lock()
|
|||
|
defer s.mutex.Unlock()
|
|||
|
cgsID := s.idGen()
|
|||
|
cgsVector := new(vector.Vector3).SetOn(0, pos[0], pos[1])
|
|||
|
s.cgs[cgsID] = &CGS{
|
|||
|
Position: cgsVector,
|
|||
|
Radius: radius,
|
|||
|
Frequency: frequency,
|
|||
|
}
|
|||
|
s.rtree.Insert(cgsVector) // Добавление в R-Tree
|
|||
|
}
|
|||
|
|
|||
|
// Adding UAV with its path and coverage radius
|
|||
|
func (s *Simulation) AddUAV(path map[int][]float64, speed []float64, radius float64) {
|
|||
|
s.mutex.Lock()
|
|||
|
defer s.mutex.Unlock()
|
|||
|
uavID := s.idGen()
|
|||
|
speedVector := &vector.Vector3{X: speed[0], Y: speed[1], Z: speed[2]}
|
|||
|
freq := 5030 + rand.Float64()*(5090-5030) // Генерация частоты в диапазоне 5030-5090 МГц
|
|||
|
|
|||
|
s.uavs[uavID] = &UAV{
|
|||
|
Path: path,
|
|||
|
Speed: speedVector,
|
|||
|
Freq: freq,
|
|||
|
Radius: radius,
|
|||
|
}
|
|||
|
|
|||
|
for _, p := range path {
|
|||
|
pos := new(vector.Vector3).SetOn(p[0], p[1], p[2])
|
|||
|
s.rtree.Insert(pos) // Добавление позиций в R-Tree
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
// Removing UAV
|
|||
|
func (s *Simulation) RemoveUAV(id int) {
|
|||
|
s.mutex.Lock()
|
|||
|
defer s.mutex.Unlock()
|
|||
|
delete(s.uavs, id)
|
|||
|
}
|
|||
|
|
|||
|
func (s *Simulation) GetPosition(uav *UAV, time int) *vector.Vector3 {
|
|||
|
// Determine known time slices
|
|||
|
keys := make([]int, 0)
|
|||
|
for k := range uav.Path {
|
|||
|
keys = append(keys, k)
|
|||
|
}
|
|||
|
|
|||
|
if _, exists := uav.Path[time]; exists {
|
|||
|
return new(vector.Vector3).SetOn(uav.Path[time][0], uav.Path[time][1], uav.Path[time][2])
|
|||
|
}
|
|||
|
|
|||
|
// If out of known range, use speed to estimate position
|
|||
|
if time < keys[0] || time > keys[len(keys)-1] {
|
|||
|
initialTime := keys[len(keys)-1]
|
|||
|
initialPosition := &vector.Vector3{
|
|||
|
X: uav.Path[initialTime][1],
|
|||
|
Y: uav.Path[initialTime][2],
|
|||
|
Z: uav.Path[initialTime][0],
|
|||
|
}
|
|||
|
deltaTime := float64(time - initialTime)
|
|||
|
speedComponent := uav.Speed.Scale(deltaTime)
|
|||
|
return initialPosition.Add(speedComponent)
|
|||
|
}
|
|||
|
|
|||
|
// Otherwise interpolate
|
|||
|
keys = append(keys, time)
|
|||
|
keys = utils.Sorted(keys)
|
|||
|
idx := utils.IndexOf(keys, time)
|
|||
|
temp := new(vector.Vector3)
|
|||
|
temp.RerpVectors(
|
|||
|
&vector.Vector3{X: uav.Path[keys[idx-1]][1], Y: uav.Path[keys[idx-1]][2], Z: uav.Path[keys[idx-1]][0]},
|
|||
|
&vector.Vector3{X: uav.Path[keys[idx+1]][1], Y: uav.Path[keys[idx+1]][2], Z: uav.Path[keys[idx+1]][0]},
|
|||
|
float64(time-keys[idx-1])/float64(keys[idx+1]-keys[idx-1]),
|
|||
|
)
|
|||
|
return temp.Copy()
|
|||
|
}
|
|||
|
|
|||
|
// Setting timeline configuration
|
|||
|
func (s *Simulation) SetTimeline(conf []int) {
|
|||
|
s.timeline = make([]int, 0)
|
|||
|
for i := conf[0]; i < conf[1]; i += conf[2] {
|
|||
|
s.timeline = append(s.timeline, i)
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
// Calculating Great Circle Distance (to account for Earth's curvature)
|
|||
|
func GreatCircleDistance(v1, v2 *vector.Vector3) float64 {
|
|||
|
lat1, lon1 := v1.Y*(math.Pi/180), v1.X*(math.Pi/180)
|
|||
|
lat2, lon2 := v2.Y*(math.Pi/180), v2.X*(math.Pi/180)
|
|||
|
dlat := lat2 - lat1
|
|||
|
dlon := lon2 - lon1
|
|||
|
a := math.Sin(dlat/2)*math.Sin(dlat/2) + math.Cos(lat1)*math.Cos(lat2)*math.Sin(dlon/2)*math.Sin(dlon/2)
|
|||
|
c := 2 * math.Atan2(math.Sqrt(a), math.Sqrt(1-a))
|
|||
|
R := 6371.0 // Радиус Земли в км
|
|||
|
return R * c
|
|||
|
}
|
|||
|
|
|||
|
// FreeSpacePathLoss Calculating Free Space Path Loss (FSPL)
|
|||
|
func FreeSpacePathLoss(distance float64, frequency float64) float64 {
|
|||
|
return 20*math.Log10(distance*1e3) + 20*math.Log10(frequency) - 20*math.Log10(3e8/frequency) + 32.44
|
|||
|
}
|
|||
|
|
|||
|
// Determining Line of Sight (LoS)
|
|||
|
// Function to represent Line of Sight (LoS) check considering terrain and Fresnel zones
|
|||
|
func (s *Simulation) HasLineOfSight(v1, v2 *vector.Vector3) bool {
|
|||
|
// Расчет расстояния и средней точки между v1 и v2
|
|||
|
distance := GreatCircleDistance(v1, v2)
|
|||
|
midpoint := &vector.Vector3{
|
|||
|
X: (v1.X + v2.X) / 2,
|
|||
|
Y: (v1.Y + v2.Y) / 2,
|
|||
|
Z: (v1.Z + v2.Z) / 2,
|
|||
|
}
|
|||
|
|
|||
|
// Проверка высоты на средней точке для оценки наличия препятствий
|
|||
|
key := fmt.Sprintf("%.6f,%.6f", midpoint.X, midpoint.Y)
|
|||
|
terrainHeight, exists := s.heights[key]
|
|||
|
if exists && midpoint.Z < terrainHeight {
|
|||
|
return false
|
|||
|
}
|
|||
|
|
|||
|
// Проверка первой зоны Френеля
|
|||
|
fresnelRadius := FresnelZoneRadius(distance, v1.Length(), v2.Length())
|
|||
|
if midpoint.Z < terrainHeight+fresnelRadius {
|
|||
|
return false
|
|||
|
}
|
|||
|
|
|||
|
return true
|
|||
|
}
|
|||
|
|
|||
|
// Calculate the radius of the first Fresnel zone
|
|||
|
func FresnelZoneRadius(distance, height1, height2 float64) float64 {
|
|||
|
wavelength := 3e8 / 5.06e9 // Средняя частота 5060 МГц
|
|||
|
return math.Sqrt((wavelength * distance) / (height1 + height2))
|
|||
|
}
|
|||
|
|
|||
|
// Calculating connectivity between UAVs and CGSs
|
|||
|
func (s *Simulation) CalculateConnectivity() {
|
|||
|
for _, uav := range s.uavs {
|
|||
|
for time := range s.timeline {
|
|||
|
uavPosition := s.GetPosition(uav, time)
|
|||
|
if uavPosition == nil {
|
|||
|
continue
|
|||
|
}
|
|||
|
for _, cgs := range s.cgs {
|
|||
|
distance := GreatCircleDistance(uavPosition, cgs.Position)
|
|||
|
if distance <= uav.Radius && distance <= cgs.Radius && s.HasLineOfSight(uavPosition, cgs.Position) {
|
|||
|
fspl := FreeSpacePathLoss(distance, uav.Freq*1e6) // Используем частоту UAV в МГц
|
|||
|
snr := 100.0 / (fspl + 1) // Примерный SNR для расчетов
|
|||
|
modulation := s.CalculateModulation(distance)
|
|||
|
dataRate := s.CalculateDataRate(modulation, snr)
|
|||
|
|
|||
|
fmt.Printf("\u0412ремя %d: UAV на позиции (%.2f, %.2f, %.2f) может подключиться к CGS на позиции (%.2f, %.2f, %.2f) с потерями %.2f дБ, модуляция: %s, частота CGS: %.2f МГц, скорость передачи: %.2f Mbps\n",
|
|||
|
time, uavPosition.X, uavPosition.Y, uavPosition.Z,
|
|||
|
cgs.Position.X, cgs.Position.Y, cgs.Position.Z, fspl, modulation, cgs.Frequency, dataRate)
|
|||
|
}
|
|||
|
}
|
|||
|
s.results = append(s.results, s.CollectResult(time))
|
|||
|
}
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
// Calculate Modulation Scheme Based on Distance
|
|||
|
func (s *Simulation) CalculateModulation(distance float64) string {
|
|||
|
if distance < 500 {
|
|||
|
return "64-QAM"
|
|||
|
} else if distance < 1000 {
|
|||
|
return "16-QAM"
|
|||
|
} else if distance < 2000 {
|
|||
|
return "QPSK"
|
|||
|
}
|
|||
|
return "BPSK"
|
|||
|
}
|
|||
|
func (s *Simulation) CalculateDataRate(modulation string, snr float64) float64 {
|
|||
|
var spectralEfficiency float64
|
|||
|
switch modulation {
|
|||
|
case "64-QAM":
|
|||
|
spectralEfficiency = 6.0
|
|||
|
case "16-QAM":
|
|||
|
spectralEfficiency = 4.0
|
|||
|
case "QPSK":
|
|||
|
spectralEfficiency = 2.0
|
|||
|
case "BPSK":
|
|||
|
spectralEfficiency = 1.0
|
|||
|
default:
|
|||
|
spectralEfficiency = 0.5
|
|||
|
}
|
|||
|
|
|||
|
// Подсчет пропускной способности (в Mbps)
|
|||
|
bandwidth := 20.0 // Пропускная способность канала в МГц
|
|||
|
dataRate := bandwidth * spectralEfficiency * math.Log2(1+snr)
|
|||
|
return dataRate
|
|||
|
}
|
|||
|
|
|||
|
// CollectResult collects the current state of the simulation at a given time
|
|||
|
func (s *Simulation) CollectResult(time int) SimulationResult {
|
|||
|
uavsResult := make(map[int]UAVResult)
|
|||
|
for id, uav := range s.uavs {
|
|||
|
position := s.GetPosition(uav, time)
|
|||
|
|
|||
|
if position != nil {
|
|||
|
uavsResult[id] = UAVResult{
|
|||
|
Position: position,
|
|||
|
Freq: uav.Freq,
|
|||
|
Radius: uav.Radius,
|
|||
|
DataRate: s.CalculateDataRate(s.CalculateModulation(GreatCircleDistance(position, s.cgs[1].Position)),
|
|||
|
100.0/(FreeSpacePathLoss(GreatCircleDistance(position, s.cgs[1].Position), uav.Freq*1e6)+1)),
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
}
|
|||
|
|
|||
|
cgsResult := make(map[int]CGSResult)
|
|||
|
for id, cgs := range s.cgs {
|
|||
|
cgsResult[id] = CGSResult{
|
|||
|
Position: cgs.Position,
|
|||
|
Radius: cgs.Radius,
|
|||
|
Frequency: cgs.Frequency,
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
return SimulationResult{
|
|||
|
Time: time,
|
|||
|
UAVs: uavsResult,
|
|||
|
CGSs: cgsResult,
|
|||
|
}
|
|||
|
}
|
|||
|
|
|||
|
// GetResultsJSON returns the results of the simulation as a JSON string
|
|||
|
func (s *Simulation) GetResultsJSON() (string, error) {
|
|||
|
jsonData, err := json.MarshalIndent(s.results, "", " ")
|
|||
|
if err != nil {
|
|||
|
return "", err
|
|||
|
}
|
|||
|
return string(jsonData), nil
|
|||
|
}
|
|||
|
|
|||
|
// Test function to create example data and run the simulation
|
|||
|
func TestSimulation() []byte {
|
|||
|
sim := NewSimulation()
|
|||
|
|
|||
|
// Add CGSs (Ground Stations)
|
|||
|
sim.AddCGS([]float64{55.7558, 37.6176}, 100, 5060)
|
|||
|
sim.AddCGS([]float64{56.7558, 37.7176}, 80, 5080)
|
|||
|
sim.AddCGS([]float64{34.0522, -118.2437}, 120, 5035)
|
|||
|
|
|||
|
// Add UAVs (drones) with paths and speed
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1000, 55.7558, 37.6176},
|
|||
|
10: {1100, 56.0, 38.0},
|
|||
|
20: {2200, 256.5, 78.5},
|
|||
|
}, []float64{10, 0.1, 0.1}, 100)
|
|||
|
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1050, 55.7558, 37.6176},
|
|||
|
10: {1100, 56.0, 38.0},
|
|||
|
20: {1200, 56.5, 38.5},
|
|||
|
}, []float64{100, 0.1, 0.1}, 50)
|
|||
|
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1010, 55.7558, 37.6176},
|
|||
|
10: {1100, 56.0, 38.0},
|
|||
|
20: {1200, 56.5, 38.5},
|
|||
|
}, []float64{10, 0.1, 0.1}, 100)
|
|||
|
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1500, 40.7128, -74.0060},
|
|||
|
10: {1600, 41.0, -73.5},
|
|||
|
20: {1700, 41.5, -73.0},
|
|||
|
}, []float64{12, 0.15, 0.1}, 100)
|
|||
|
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1400, 34.0522, -118.2437},
|
|||
|
10: {1500, 34.5, -118.0},
|
|||
|
20: {1600, 35.0, -117.5},
|
|||
|
}, []float64{15, 0.2, 0.1}, 100)
|
|||
|
|
|||
|
sim.AddUAV(map[int][]float64{
|
|||
|
0: {1300, 48.8566, 2.3522},
|
|||
|
10: {1350, 49.0, 2.5},
|
|||
|
20: {1400, 49.5, 3.0},
|
|||
|
}, []float64{8, 0.1, 0.05}, 100)
|
|||
|
|
|||
|
// Set timeline configuration [start, end, step]
|
|||
|
sim.SetTimeline([]int{0, 50, 1})
|
|||
|
|
|||
|
// Run simulation and print results
|
|||
|
fmt.Println("Результаты симуляции:")
|
|||
|
for _, t := range sim.timeline {
|
|||
|
fmt.Printf("\nВремя %d секунд:\n", t)
|
|||
|
for uavID, uav := range sim.uavs {
|
|||
|
position := sim.GetPosition(uav, t)
|
|||
|
if position != nil {
|
|||
|
fmt.Printf("UAV %d на позиции (%.2f, %.2f, %.2f)\n", uavID, position.X, position.Y, position.Z)
|
|||
|
}
|
|||
|
}
|
|||
|
sim.CalculateConnectivity()
|
|||
|
}
|
|||
|
result, err := json.Marshal(sim.results)
|
|||
|
if err != nil {
|
|||
|
panic(err)
|
|||
|
}
|
|||
|
return result
|
|||
|
}
|