wgs84转ENU的matlab代码

以下是wgs84转ENU的matlab代码：

function [east,north,up] = wgs84_to_enu(lat,lon,height,lat_ref,lon_ref,height_ref)

% Convert WGS84 coordinates to ENU coordinates
% Inputs:
% lat, lon, height: WGS84 coordinates of the point of interest (in degrees and meters)
% lat_ref, lon_ref, height_ref: WGS84 coordinates of the reference point (in degrees and meters)
% Outputs:
% east, north, up: ENU coordinates of the point of interest (in meters)

% Convert degrees to radians
lat = deg2rad(lat);
lon = deg2rad(lon);
lat_ref = deg2rad(lat_ref);
lon_ref = deg2rad(lon_ref);

% WGS84 constants
a = 6378137; % semi-major axis (m)
f = 1/298.257223563; % flattening
b = a*(1-f); % semi-minor axis (m)
e = sqrt((a^2-b^2)/a^2); % eccentricity

% Calculate the radius of curvature in the prime vertical
N = a/sqrt(1-e^2*sin(lat_ref)^2);

% Calculate the ENU coordinates
east = (sin(lon-lon_ref)*cos(lat)*N + cos(lon-lon_ref)*cos(lat)*N + sin(lat-lat_ref)*height)*1000;
north = (-sin(lon-lon_ref)*sin(lat)*N + cos(lon-lon_ref)*sin(lat)*N + cos(lat-lat_ref)*height)*1000;
up = (cos(lon-lon_ref)*cos(lat)*height + sin(lon-lon_ref)*cos(lat)*height + sin(lat-lat_ref)*N)*1000;

end

请注意，此代码仅适用于计算地球表面上的点的ENU坐标。如果您需要处理大范围的距离或高程变化，则需要使用更复杂的模型和算法，例如大地坐标系（Geodetic Coordinate System）和高斯-克吕格投影（Gauss-Kruger Projection）。