matrix calculation reshaping and diferences between datum

Question

0 个投票

Coordinate_Systems.m

i dont know how to transform the Z dimension into a 2x2 matriz at least i dont understand how to make it run

7 个评论
显示 5更早的评论隐藏 5更早的评论

Voss 2024-1-2

Where would the country border come from?

I don't know what's going on either; you'll need to step through your code and check that each variable seems reasonable (correct size, reasonable values) as it runs, and find where it goes wrong and fix it.

Vasco 2024-1-2

i checked the values, the values of the other variables check, i just dont get what is happening maybe is it because the data doesnt come sorted?

the figure im making is an attempt to answer this question :"Assess the differences between ETRS89 and the local datums (D73 and DLx) regarding geodesic coordinates (Δφ, Δλ, Δh)."

请先登录，再进行评论。

请先登录，再回答此问题。

Follow Question

Answer 1

Cris LaPierre 2024-1-2

在 MATLAB Online 中打开

0 个投票

RGN_ETRS89_geod.txt

I'd use scatteredInterpolant probably.

Also, you can simplify your code a little by accessing the data in your tables more directly. See here: https://www.mathworks.com/help/matlab/matlab_prog/access-data-in-a-table.html

Here's the relevant pieces of code. Note that I am not claiming to solve your assignment for you. Just get you started.

%Inserção dos dados e criação de variáveis

file_in=readtable("RGN_ETRS89_geod.txt");

lat_dec = dms2degrees(file_in{:,["Latd","Latm","Lats"]});

lon_dec = -1*dms2degrees(file_in{:,["Lond","Lonm","Lons"]});

alt_orto = file_in.h;

alt_elips = file_in.H;

% scatterm(lat_dec,lon_dec)

k = boundary(lat_dec,lon_dec);

geoshow(lat_dec(k),lon_dec(k))

%Elipsoide de referência para ETRS89: GRS80

a=6378137;

f=1/298.257222101;

b=a*(1-f);

e=(sqrt(a^2-b^2))/a;

N=a./sqrt(1-e^2.*(sind(lat_dec)).^2);

%Cálculo de X, Y e Z

X=(N+alt_elips).*cosd(lat_dec).*cosd(lon_dec);

Y=(N+alt_elips).*cosd(lat_dec).*sind(lon_dec);

Z=(N.*(1-e^2)+alt_elips).*sind(lat_dec);

%Exercício 2

%2.1

Tx=230.994;

Ty=-102.591;

Tz=-25.199;

omegax=deg2rad(-0.633/3600);

omegay=deg2rad(0.239/3600);

omegaz=deg2rad(-0.900/3600);

D=-1.950*10^-6;

matriz=[1 -omegaz omegay; omegaz 1 -omegax; -omegay omegax 1];

for i=1:1668

D73=[Tx;Ty;Tz]+(D+1)*matriz*[X(i);Y(i);Z(i)];

XD73(i,1)=D73(1);

YD73(i,1)=D73(2);

ZD73(i,1)=D73(3);

end

%2.2

Tx=283.088;

Ty=70.693;

Tz=-117.445;

omegax=deg2rad(1.157/3600);

omegay=deg2rad(-0.059/3600);

omegaz=deg2rad(0.652/3600);

D=4.058*10^-6;

matriz=[1 -omegaz omegay; omegaz 1 -omegax; -omegay omegax 1];

for i=1:1668

DLx=[Tx;Ty;Tz]+(D+1)*matriz*[X(i);Y(i);Z(i)];

XDLx(i,1)=DLx(1);

YDLx(i,1)=DLx(2);

ZDLx(i,1)=DLx(3);

end

%Elipsóide de referência para o Datum 73 e Datum Lisboa - Hayford

a_D73= 6378388;

f_D73=1/297;

b_D73=a_D73*(1-f_D73);

e_D73=(sqrt(a_D73^2-b_D73^2))/a_D73;

%Datum Lisboa

%Determinação da latitude

%Estimativa inicial

p=sqrt(XDLx.^2+YDLx.^2);

lat_0=atand(ZDLx./p);

N_0=a_D73./sqrt(1-e_D73^2.*(sind(lat_0)).^2);

%Primeira iteração

lat_1=atand((ZDLx+e_D73^2.*N_0.*sind(lat_0))./p);

N_1=a_D73./sqrt(1-e_D73^2.*(sind(lat_1)).^2);

%Segunda iteração

lat_2=atand((ZDLx+e_D73^2.*N_1.*sind(lat_1))./p);

N_2=a_D73./sqrt(1-e_D73^2.*(sind(lat_2)).^2);

%Terceira iteração

lat_3_DLx=atand((ZDLx+e_D73^2.*N_2.*sind(lat_2))./p);

%DLX figura e D73 figura

delta_lat_Dlx = lat_dec - lat_3_DLx;

% Create a function for predicting delta_lat_Dlx based on lat & lon data

F = scatteredInterpolant(lon_dec,lat_dec,delta_lat_Dlx);

% Create a 2d grid of Lat and Lon values

[LON,LAT] = meshgrid(-10:0.1:-6,36:0.05:43);

% Predict gridded values for delta_lat_Dlx

Z = F(LON,LAT);

% set values outside portugal to nan (so they do not appear in contour plot

in = inpolygon(LAT,LON,lat_dec(k),lon_dec(k));

Z(~in) = nan;

% Create a contour plot

hold on

contourm(LAT,LON,Z,10,'r','ShowText','on')

hold off

6 个评论
显示 4更早的评论隐藏 4更早的评论

Vasco 2024-1-3

编辑：Cris LaPierre 2024-1-3

在 MATLAB Online 中打开

First of all thank you for your suport.

I tried doing it for the lat, lon and h

the figure appears but the limits of the country dont Im gonna send the full code again if you could help me fix the problem, it would be live saving.

close all;

clear all;

clc;

%Exercício 1

%Inserção dos dados e criação de variáveis

file_in=readtable("RGN_ETRS89_geod.txt");

lat_dec = dms2degrees(file_in{:,["Latd","Latm","Lats"]});

lon_dec = -1*dms2degrees(file_in{:,["Lond","Lonm","Lons"]});

alt_orto = file_in.h;

alt_elips = file_in.H;

%Transformação dos dados em decimais

k = boundary(lat_dec,lon_dec);

geoshow(lat_dec(k),lon_dec(k));

%Elipsoide de referência para ETRS89: GRS80

a=6378137;

f=1/298.257222101;

b=a*(1-f);

e=(sqrt(a^2-b^2))/a;

%Cálculo de N (grande normal)

N=a./sqrt(1-e^2.*(sind(lat_dec)).^2);

%Cálculo de X, Y e Z

X=(N+alt_elips).*cosd(lat_dec).*cosd(lon_dec);

Y=(N+alt_elips).*cosd(lat_dec).*sind(lon_dec);

Z=(N.*(1-e^2)+alt_elips).*sind(lat_dec);

%%

%Exercício 2

%2.1

Tx=230.994;

Ty=-102.591;

Tz=-25.199;

omegax=deg2rad(-0.633/3600);

omegay=deg2rad(0.239/3600);

omegaz=deg2rad(-0.900/3600);

D=-1.950*10^-6;

matriz=[1 -omegaz omegay; omegaz 1 -omegax; -omegay omegax 1];

for i=1:1668

D73=[Tx;Ty;Tz]+(D+1)*matriz*[X(i);Y(i);Z(i)];

XD73(i,1)=D73(1);

YD73(i,1)=D73(2);

ZD73(i,1)=D73(3);

end

%%

%2.2

Tx=283.088;

Ty=70.693;

Tz=-117.445;

omegax=deg2rad(1.157/3600);

omegay=deg2rad(-0.059/3600);

omegaz=deg2rad(0.652/3600);

D=4.058*10^-6;

matriz=[1 -omegaz omegay; omegaz 1 -omegax; -omegay omegax 1];

for i=1:1668

DLx=[Tx;Ty;Tz]+(D+1)*matriz*[X(i);Y(i);Z(i)];

XDLx(i,1)=DLx(1);

YDLx(i,1)=DLx(2);

ZDLx(i,1)=DLx(3);

end

%%

%Elipsóide de referência para o Datum 73 e Datum Lisboa - Hayford

a_D73= 6378388;

f_D73=1/297;

b_D73=a_D73*(1-f_D73);

e_D73=(sqrt(a_D73^2-b_D73^2))/a_D73;

%Cálculo da longitude

lon_D73=atand(YD73(:,1)./XD73(:,1));

%Determinação da latitude

%Estimativa inicial

p=sqrt((XD73(:,1)).^2+(YD73(:,1)).^2);

lat_0_73=atand(ZD73(:,1)./p);

N_0_73=a_D73./sqrt(1-e_D73^2.*(sind(lat_0_73)).^2);

%Primeira iteração

lat_1_73=atand((ZD73(:,1)+e_D73^2.*N_0_73.*sind(lat_0_73))./p);

N_1_73=a_D73./sqrt(1-e_D73^2.*(sind(lat_1_73)).^2);

%Segunda iteração

lat_2_73=atand((ZD73(:,1)+e_D73^2.*N_1_73.*sind(lat_1_73))./p);

N_2_73=a_D73./sqrt(1-e_D73^2.*(sind(lat_2_73)).^2);

%Terceira iteração

lat_3_D73=atand((ZD73(:,1)+e_D73^2.*N_2_73.*sind(lat_2_73))./p);

N_3_73=a_D73./sqrt(1-e_D73^2.*(sind(lat_3_D73)).^2);

%Os valores já não se alteram da segunda para a terceira iteração, fazendo

%desta a que fornece o valor definitivo da latitude.

%Cálculo de h

h_D73=(p./cosd(lat_3_D73))-N_3_73;

%%

%Datum Lisboa

%Cálculo da longitude

lon_DLx=atand(YDLx./XDLx);

%Determinação da latitude

%Estimativa inicial

p=sqrt(XDLx.^2+YDLx.^2);

lat_0_lx=atand(ZDLx./p);

N_0_lx=a_D73./sqrt(1-e_D73^2.*(sind(lat_0_lx)).^2);

%Primeira iteração

lat_1_lx=atand((ZDLx+e_D73^2.*N_0_lx.*sind(lat_0_lx))./p);

N_1_lx=a_D73./sqrt(1-e_D73^2.*(sind(lat_1_lx)).^2);

%Segunda iteração

lat_2_lx=atand((ZDLx+e_D73^2.*N_1_lx.*sind(lat_1_lx))./p);

N_2_lx=a_D73./sqrt(1-e_D73^2.*(sind(lat_2_lx)).^2);

%Terceira iteração

lat_3_DLx=atand((ZDLx+e_D73^2.*N_2_lx.*sind(lat_2_lx))./p);

N_3_lx=a_D73./sqrt(1-e_D73^2.*(sind(lat_3_DLx)).^2);

%Os valores já não se alteram da segunda para a terceira iteração, fazendo

%desta a que fornece o valor definitivo da latitude.

%Cálculo de h

h_DLx=(p./cosd(lat_3_DLx))-N_3_lx;

%%

%----------------

%----- etrs89 and datum lisboa

%DLX figura e D73 figura

delta_lat_Dlx = (lat_dec - lat_3_DLx).*3600;

delta_lon_Dlx = (lon_dec - lon_D73).*3600;

delta_h_Dlx = (alt_elips-h_D73);

%%

figure(1);

% Create a function for predicting delta_lat_Dlx based on lat & lon data

F = scatteredInterpolant(lon_dec,lat_dec,delta_lat_Dlx);

% Create a 2d grid of Lat and Lon values

[LON,LAT] = meshgrid(-10:0.1:-6,36:0.05:43);

% Predict gridded values for delta_lat_Dlx

Z = F(LON,LAT);

% set values outside portugal to nan (so they do not appear in contour plot

in = inpolygon(LAT,LON,lat_dec(k),lon_dec(k));

Z(~in) = nan;

% Create a contour plot

hold on

contourm(LAT,LON,Z,10,'r','ShowText','on')

hold off

%%

figure(2);

% Create a function for predicting delta_lat_Dlx based on lat & lon data

G = scatteredInterpolant(lon_dec,lat_dec,delta_lon_Dlx);

% Create a 2d grid of Lat and Lon values

[LON_1,LAT_1] = meshgrid(-10:0.1:-6,36:0.05:43);

% Predict gridded values for delta_lat_Dlx

V = G(LON_1,LAT_1);

% set values outside portugal to nan (so they do not appear in contour plot

in_1 = inpolygon(LAT_1,LON_1,lat_dec(k),lon_dec(k));

V(~in_1) = nan;

% Create a contour plot

hold on

contourm(LAT_1,LON_1,V,10,'r','ShowText','on')

hold off

%%

figure(3);

% Create a function for predicting delta_lat_Dlx based on lat & lon data

L = scatteredInterpolant(lon_dec,lat_dec,delta_h_Dlx);

% Create a 2d grid of Lat and Lon values

[LON_2,LAT_2] = meshgrid(-10:0.1:-6,36:0.05:43);

% Predict gridded values for delta_lat_Dlx

O = L(LON_2,LAT_2);

% set values outside portugal to nan (so they do not appear in contour plot

in = inpolygon(LAT_2,LON_2,lat_dec(k),lon_dec(k));

O(~in_2) = nan;

Unrecognized function or variable 'in_2'.

% Create a contour plot

hold on

contourm(LAT_2,LON_2,O,10,'r','ShowText','on')

hold off

Vasco 2024-1-3

thank you very much its fixed

i try to do it with subplot but theres some type of error because it loads the lines but it doesnt load the values nor the borders.

Cris LaPierre 2024-1-4

I wouldn't expect subplot to fix the issue. Please try what I suggested.

请先登录，再进行评论。

matrix calculation reshaping and diferences between datum

7 个评论
显示 5更早的评论隐藏 5更早的评论

采纳的回答

6 个评论
显示 4更早的评论隐藏 4更早的评论

更多回答（0 个）

类别

标签

Community Treasure Hunt

matrix calculation reshaping and diferences between datum

7 个评论 显示 5更早的评论 隐藏 5更早的评论

采纳的回答

6 个评论 显示 4更早的评论 隐藏 4更早的评论

更多回答（0 个）

类别

标签

另请参阅

Community Treasure Hunt

7 个评论
显示 5更早的评论隐藏 5更早的评论

6 个评论
显示 4更早的评论隐藏 4更早的评论