This commit is contained in:
shuibing811
2025-01-24 16:11:34 +08:00
parent 07c9056cc4
commit d0932b7d18
2 changed files with 190 additions and 13 deletions

View File

@@ -0,0 +1,125 @@
import marimo
__generated_with = "0.10.16"
app = marimo.App(width="medium")
@app.cell
def _():
import numpy as np
import matplotlib.pyplot as plt
import marimo as mo
import define
return define, mo, np, plt
@app.cell
def _(lla2ecef, np):
class Station:
def __init__(self, lat, lon, alt):
self.lla = np.array([lat, lon, alt])
self.xyz = np.array(lla2ecef(*self.lla))
class Target:
def __init__(self, lat, lon, alt):
self.lla = np.array([lat, lon, alt])
self.xyz = np.array(lla2ecef(*self.lla))
class Radar:
def __init__(self, lat, lon, alt, /, *, pri_info, rf_info, pt, gt, pw):
self.lla = np.array([lat, lon, alt])
self.xyz = np.array(lla2ecef(*self.lla))
self.pri_info = pri_info
self.rf_info = rf_info
self.pt = pt
self.gt = gt
self.pw = pw
class RadarPdw:
def __init__(self):
pass
return Radar, RadarPdw, Station, Target
@app.cell
def _(Station, Target):
stations = [Station(31.2304, 121.4737, 0), Station(31.5200, 121.6400, 0)]
targets = [
Target(31.2550, 121.4890, 0), # 目标1 陆
Target(31.2800, 121.5200, 0), # 目标2 陆
Target(31.1900, 121.4800, 0), # 目标3 陆
Target(31.2700, 121.5500, 0), # 目标4 陆
# Target(31.2000, 121.6000, 0), # 目标5 陆
Target(30.9500, 121.6500, 0), # 目标6 海
Target(30.8500, 121.7200, 0), # 目标7 海
Target(30.7800, 121.8000, 0), # 目标8 海
Target(30.9200, 121.8500, 0), # 目标9 海
Target(30.6700, 121.9000, 0), # 目标10 海
]
return stations, targets
@app.cell
def _():
return
@app.cell
def _():
import pyproj
def lla2ecef(lat,lon,alt):
"""
将经纬度高度LLA转换为地心地固坐标系ECEF坐标。
参数:
lon (float): 经度,单位为度。
lat (float): 纬度,单位为度。
alt (float): 海拔,单位为米。
返回:
tuple: 包含 ECEF 坐标 (x, y, z) 的元组,单位为米。
"""
# 定义源坐标系统WGS84 经纬度)
lla = pyproj.CRS("EPSG:4326")
# 定义目标坐标系统ECEF
ecef = pyproj.CRS("EPSG:4978")
# 创建转换器对象
transformer = pyproj.Transformer.from_crs(lla, ecef)
# 进行坐标转换
x, y, z = transformer.transform(lat, lon, alt, radians=False)
return x, y, z
def lla2enu(lat, lon, alt, lat0, lon0, alt0):
"""
将经度、纬度和高度LLA坐标转换为东 - 北 - 天ENU坐标。
参数:
lon (float): 待转换点的经度(度)
lat (float): 待转换点的纬度(度)
alt (float): 待转换点的高度(米)
lon0 (float): 参考点的经度(度)
lat0 (float): 参考点的纬度(度)
alt0 (float): 参考点的高度(米)
返回:
tuple: 包含 E、N、U坐标的元组
"""
# 定义 WGS84 经纬度坐标系统
lla = pyproj.CRS("EPSG:4326")
# 定义局部的 ENU 坐标系统
enu = pyproj.CRS.from_string(f"+proj=aeqd +lat_0={lat0} +lon_0={lon0} +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs")
# 创建坐标转换器
transformer = pyproj.Transformer.from_crs(lla, enu)
# 执行坐标转换
e, n, u = transformer.transform(lon, lat, alt)
u = u - alt0
return e, n, u
return lla2ecef, lla2enu, pyproj
if __name__ == "__main__":
app.run()

78
main.py
View File

@@ -7,11 +7,10 @@ app = marimo.App(width="medium")
@app.cell
def _():
import numpy as np
import random
import math
import matplotlib.pyplot as plt
import marimo as mo
return math, mo, np, plt, random
return math, mo, np, plt
@app.cell
@@ -157,7 +156,6 @@ def _(Radar, Station, Target):
),
]
assert(len(targets) == len(radars))
return radars, stations, targets
@@ -209,16 +207,56 @@ def _(
return (generate_pdw,)
app._unparsable_cell(
r"""
def locate():
stations_enu = np.empty((len(stations),3))
for i in len(stations):
stations_enu[i,:] = np.array(lla2enu(stations[i].lla,stations[0]))
for i in len(radar_pdws):
doa = radar_pdws[i].doa
erro = np.empty(len(doa))
location_estimated = np.empty(len(doa))
for j in len(doa):
ned = least_square(doa[j],stations_enu(:,0),stations_enu(:,1))
radi = ned2lla([ned[0],ned[1],0],stations[0].lla)
erro[j] = distance(targets[0].lla[0],targets[0].lla[1],radi[0],radi[1])
location_estimated[j] = np.array([radi(0),radi(1),0])
for j in len(radar_pdws[i].pri):
radar_pdws[i].location[j,:] = location_estimated[j,:]
radar_pdws[i].erro[j,:] = erro[j]
return radar_pdws
""",
name="_"
)
@app.cell
def _(generate_pdw):
def _(generate_pdw, locate, 保存radar_pdws):
radar_pdws = generate_pdw()
return (radar_pdws,)
#保存radar_pdws
radar_pdws_located = locate()
保存radar_pdws
return radar_pdws, radar_pdws_located
@app.cell
def _(math, np, random):
def _(math, np):
def distance(lat1,lon1,lat2,lon2):
pass
def least_square(doa,e,n):
doa_rad = np.deg2rad(doa)
l = len(doa_rad)
ee = e[0:l]
nn = n[0:l]
fmtx = ee - nn * np.tan(doa_rad)
hmtx = np.hstack(-np.tan(doa_rad).T,np.ones(l).T)
ned = (hmtx.T @ hmtx) ** -1 * hmtx.T * fmtx.T
return ned
def calc_propagation_time(radar_pos,station_pos):
return np.linalg.norm(radar_pos-station_pos)/3e8
@@ -241,7 +279,7 @@ def _(math, np, random):
elif rf_info["type"] == "Agility": #捷变
rf = np.empty(num)
for i in range(num):
rf[i] = rf_info["center"] + rf_info["bw"] * (2 * random.random() - 1)
rf[i] = rf_info["center"] + rf_info["bw"] * (2 * np.random.rand() - 1)
return rf
else:
raise ValueError("未知的rf类型")
@@ -260,7 +298,7 @@ def _(math, np, random):
m = 10
pri = np.full(m,0)
for i in range(m):
gamma = (0.01 + (0.2 - 0.01) * random.random()) * [-1.0,1.0][np.random.randint(2)]
gamma = (0.01 + (0.2 - 0.01) * np.random.rand()) * [-1.0,1.0][np.random.randint(2)]
pri[i] = pri_value + gamma * pri_value
pri = np.tile(pri, (1, math.ceil(n/m)))
pri = pri[0:n]
@@ -296,7 +334,14 @@ def _(math, np, random):
return toa,pri
else:
raise ValueError("未知的pri类型")
return calc_doa_north, calc_propagation_time, generate_pri, generate_rf
return (
calc_doa_north,
calc_propagation_time,
distance,
generate_pri,
generate_rf,
least_square,
)
@app.cell
@@ -354,7 +399,9 @@ def _():
u = u - alt0
return e, n, u
return lla2ecef, lla2enu, pyproj
def ned2lla(ned,lla0):
pass
return lla2ecef, lla2enu, ned2lla, pyproj
@app.cell
@@ -373,8 +420,13 @@ def _(lla2enu):
@app.cell
def _():
return
def _(np):
a = np.array([1,2],ndmin = 2).transpose()
b = np.array([3,4],ndmin = 2).transpose()
np.hstack((a,b))
np.empty((2,3))
return a, b
if __name__ == "__main__":