基于非线性最小二乘法 (NLLS) 分析寻找给定伪距数据的无人机轨迹附matlab代码

✅作者简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，matlab项目合作可私信。

🍎个人主页：Matlab科研工作室

🍊个人信条：格物致知。

更多Matlab仿真内容点击👇

智能优化算法       神经网络预测       雷达通信       无线传感器        电力系统

信号处理              图像处理               路径规划       元胞自动机        无人机 

⛄ 内容介绍

这个项目由两个主要任务组成。第一种是使用非线性最小二乘法（NLLS）分析来确定给定伪距数据的无人机轨迹，并通过将其与真实无人机轨迹匹配来验证这一结果（也给出了）。第二个是计算计算的无人机轨迹的精度稀释（DOP）。

DOP是对任何给定GPS位置数据的精度的测量。五种不同类型的

将计算1个DOP：几何DOP、位置DOP、水平DOP、垂直DOP和时间DOP。它们是x、y、z、t（时钟偏置）分量Vx、Vy、Vz和Vt的方差的函数。

几何DOP是由于卫星在太空中的几何方向而产生的误差。位置DOP是无人机三维位置的误差。水平DOP是指无人机水平位置的误差，而垂直DOP是无人机高度的误差。时间DOP是由于时间精度造成的误差。

DOP值越低，无人机的瞬时位置测量就越好。良好的GDOP或PDOP被认为＜5，而PDOP的值＞10则非常差。HDOP预计在2-3之间。由于所有卫星都在接收器上方，预计VDOP将高于HDOP，而对于水平坐标，数据是从无人机的所有侧面接收的。

GPS定位问题

假设接收器位于二维地面上。两颗全球定位系统卫星可以用来精确定位其在地球上的位置，前提是如果这两颗卫星的范围由一个圆圈表示，那么将有两个交点，其中只有一个位于地球上。这没有给出范围误差。假设GPS GPS卫星存在距离误差，两个距离圆的交点现在将发生偏移，导致接收器在地球上的位置估计不正确。

这个问题可以通过让第三个GPS卫星解决接收器的位置来解决。假设所有GPS卫星具有相同的时钟/距离误差，则所有三个距离圆都将在接收器位置相交。这可以用来计算时钟偏置，从而计算范围误差。类似地，在三维空间中，我们需要四颗卫星来精确定位接收器的位置。

方法论

这项任务使用了两个给定的数据文件。第一个是GPS伪距F1.txt。该文件包含从无人机到所有可见卫星计算的伪距。第二个是UAVPostion F1.txt。该文件包含无人机的真实位置。图2中的流程图描述了根据给定伪距数据确定无人机轨迹的过程。该过程使用了非线性最小二乘回归。

定义f（x0）和H的方程式可在附录A部分找到。上图中的流程图。向量X0的形式为X0=[x，y，z，cbU]，其中cbU是时钟偏置。后者是测量地面GPS接收器中设置的时间。

理想情况下，如果地球上的GPS接收器和GPS卫星都有同步的原子钟，那么只需要三颗卫星就可以精确定位接收器的位置。然而，一个原子钟的价格高达20万美元，这使得它过于昂贵，无法在民用市场上广泛使用。或者，可以用一个便宜的接收器进行简单的权衡，这一次需要四颗卫星来确定接收器的位置。额外的卫星用于获得时钟偏移，即用于计算偏移时间以校正接收器时钟的偏移距离。

⛄ 部分代码

% Clear workspace, close all open figure & clear the command window

clear all

close all;

clc;

% Extract the data from the text file

% Add other folders to path

addpath('../data', '../lib/', '../lib/conversion');

% Load constants

constants();

% Initialise Ground Station Position

% Defining Sydney Ground Station Coordinates

lat        = deg2rad(-34.76);

long       = deg2rad(150.03);

alt        = 680;

pos_llh_gs = [lat;long;alt];

% File Name

uav_data_fpath = 'UAVPosition_F1.txt';

% Import pseudorange data

pseudo_data = importdata('GPS_pseudorange_F1.txt');

% Load ECEF position values of satellites

load ECEFPos

%% Categorising time values

% Vernal Equinox time

equinox_time = 7347737.336;

% Store time data

times     = pseudo_data(:,1) - equinox_time;

% Obtain the unique time values

% nTimes is an array containing the total number of occurences of single time

% value

% timeVal is an array containing all unique timevalues

[total_occurences time_values] = hist(times(:),unique(times));

% Cummulative Time Array

cummulative_times = cumsum(total_occurences);

% Increment time values

% timeVal = timeVal + 1;

%% Observed UAV Positions

% Obtain the polar coordinates of UAV w.r.t Ground Station

[pos_UAV_Cart_obs,pos_UAV_Pol_Obs, pos_ECEF_gs, UAVVel, UAVPos] = UAVpolarCoor(time_values,...

                                         cummulative_times, pseudo_data, ECEFPos, pos_llh_gs);

%% Extract true UAV Position data from text files

% % Obtain polar coordinates of true measurements

[pos_UAV_Pol_True, pos_UAV_Cart_True] = extractUAVtrue(uav_data_fpath,equinox_time);

%% Catalog the position of the satellites during UAV Tracking

% Find satellite position w.r.t ground station

[pos_Sat_Pol_Obs] = findSatPos(ECEFPos, pos_ECEF_gs, pos_llh_gs);

%% Calculate the DOP for each time value

% Calculate the best and worst satellite configurations according to DOP

[bestSatConfig, worstSatConfig, DOP, maxDOPIndex, minDOPIndex] = findDOP(UAVVel, pos_Sat_Pol_Obs);

%% Save relevant data as a structured array for plotting

plotData.timeArray         = times;

plotData.pos_UAV_Cart_obs  = pos_UAV_Cart_obs;

plotData.pos_UAV_Pol_Obs   = pos_UAV_Pol_Obs;

plotData.pos_ECEF_gs       = pos_ECEF_gs;

plotData.UAVVel            = UAVVel;

plotData.UAVPos            = UAVPos;

plotData.pos_UAV_Pol_True  = pos_UAV_Pol_True;

plotData.pos_UAV_Cart_True = pos_UAV_Cart_True;

plotData.pos_Sat_Pol_Obs   = pos_Sat_Pol_Obs;

plotData.bestSatConfig     = bestSatConfig;

plotData.worstSatConfig    = worstSatConfig;

plotData.DOP               = DOP;

plotData.maxDOPIndex       = maxDOPIndex;

plotData.minDOPIndex       = minDOPIndex;

plotData.nTimes            = total_occurences;

plotData.pseudo_data       = pseudo_data;

plotData.cumArray          = cummulative_times;

%% Plot graphs

plotq1B(plotData)

⛄ 运行结果

⛄ 参考文献

⛳️ 代码获取关注我

❤️部分理论引用网络文献，若有侵权联系博主删除

❤️ 关注我领取海量matlab电子书和数学建模资料