How to properly create for loop with if statement
이전 댓글 표시
Hi, how do you create a for or while loop with an if statement which adds 360 to an existing variable l, if l is negative?
댓글 수: 1
The simple and efficient approach:
L = -180:20:180
L = mod(L,360)
채택된 답변
It is necessary to use the index on ‘L’ on both sides of the equation.
L = -180:20:180
for l = 1:length(L)
if L(l) < 0
L(l) = L(l) + 360;
end
end
L
A one-line version (avoiding the loop) would be —
L = -180:20:180
L(L<0) = L(L<0) + 360
The one-line version is more efficient.
.
댓글 수: 91
Star Strider
2023년 4월 2일
As always, my pleasure!
Light
2023년 4월 2일
Star Strider
2023년 4월 2일
I am not certain where it would be applied, however if is appropriate to the angle, then yes (with appropriate changes to the variables if necessary).
If you have the Mapping Toolbox, there are several functions, such as wrapTo360 and others that could be applicable here.
‘Does wrapTo360 convert the negative longitudes to positive ones?’
It’s easiest to do the experiment since it’s available here —
L = -180:20:180;
Lw = wrapTo360(L)
So, yes.
‘... but I meant, since L was already declared as atan2d(Y, X) could it also be declared as L = -180: 20: 180 ?’
If you want to, yes. I would leave it as calculated, however, rather than re-defining it.
.
Star Strider
2023년 4월 2일
편집: Star Strider
2023년 4월 2일
You would need to establish the context there first.
I am not certain what you are asking here. I would not substitute the hard-coded ‘L’ vector for the calculated value, since the atan2 arguments could change. I would just leave it as it is.
EDIT — (2 Apr 2023 at 19:23)
I was not aware what ‘L’ is being used for. If you are using it to plot, it may be necessary to edit any other vectors being plotted as functions of it (if any are) as well. If it is a function of another variable, it could be changed without problems.
Light
2023년 4월 2일
Star Strider
2023년 4월 2일
That sounds reasonable.
What problems are you encountering?
What are the minimum and maximum angles currently being calculated? For example, would it make sense to add 180 to all of them?
Light
2023년 4월 2일
Star Strider
2023년 4월 2일
O.K.
If you are supposed to add 360 to the negative values of ‘L’ then the one-line version should work. I would put it just after ‘L’ is calculated in the atan2 call.
Light
2023년 4월 2일
Star Strider
2023년 4월 2일
I created this ‘L’ vector to test and demonstrate the code:
L = -180:20:180
That is the only reason it exists.
In your code, use the ‘L’ vector created by the atan2 call.
Star Strider
2023년 4월 2일
It would be:
L = atan2d(Y,X);
L(L<0) = L(L<0) + 360;
If I understand correctly what you want to do.
Light
2023년 4월 2일
Light
2023년 4월 2일
Star Strider
2023년 4월 2일
It would, however I only have images of it thus far. I would need the full text, ideally copied and pasted here rather than attached as .m file that I would then have to download and open.
This appears to be a homework assignment, so I will only comment on what I consider to be errors, or suggestions on improvements, and note them. I will not correct the code. If you want me to run it, I would need any data files as well.
Light
2023년 4월 2일
편집: Walter Roberson
2023년 4월 12일
Running it —
% l= longitude, b= latitude
l1 = -96.81
b1 = 32.78
% translating polar coordinates of the first city (Dallas) to cartesian coordinates
x1 = cosd (-96.81)*cosd (32.78)
y1 = sind (-96.81)*cosd (32.78)
z1 = sind (32.78)
l2 = 14.42
b2 = 50.07
% translating polar coordinates of the second city (Prague) to cartesian coordinates
x2 = cosd (14.42)*cosd (50.07)
y2 = sind (14.42)*cosd (50.07)
z2 = sind (50.07)
% calculating angular distance between the two cities (Dallas and Prague)
psi = acosd ((x1*x2)+(y1*y2)+(z1*z2))
x3 = ((x2 - x1)*cosd(psi))/ sind(psi)
y3 = ((y2 - y1)*cosd(psi))/ sind(psi)
z3 = ((z2 - z1)*cosd(psi))/ sind(psi)
% calculating the distance between the two cities (Dallas and Prague)
r = 6371;
Distance = 2*pi*r*(psi/360)
% Plotting the great circle path between the two cities
phi = linspace(0, psi, 100)
X = x1*cosd(phi) + x3*sind(phi);
Y = y1*cosd(phi) + y3*sind(phi);
Z = z1*cosd(phi) + z3*sind(phi);
% converting cartesian coordinates x, y, z to polar coordinates l, b % (longitude and latitude)
b = asind(Z)
l = atan2d (Y,X)
l(l < 0) = l(l < 0) + 360
load('topo.mat', 'topo', 'topomap1');
% load earth topography
whos topo topomap1;
% set coordinates of points
Longitude = 0:15:360; Latitude = -72:6:72;
figure
% set background parameters
contour(0:359, -89:90, topo, [0,0], 'w');
axis equal
box on
set(gca, 'Color', [0 0 1], 'XLim', [0 360], 'YLim', [-90 90], 'XTick',...
[0: 180: 360], 'YTick', [-90: 90: 90]);
%set(gca, 'Color', [0 0 1], 'XLim', [0 360], 'YLim', [-90 90]);
grid; hold on
plot(l, b, 'y.', 'linewidth', 1);
hold off
I had to insert some linefeeds to clarify the code lines, however it runs without error. I defer to you to determine if it produces the correct result.
.
Light
2023년 4월 3일
Star Strider
2023년 4월 3일
I have no idea how to assess that (although I thought the waypoint in the southern Indian Ocean was interesting).
Light
2023년 4월 3일
Star Strider
2023년 4월 3일
Takling a closer look, it would seem that perhaps half of the flightpath (for lack of a better term for it) is the mirror image of the other half, with respect to both latitude and lonogitude.
I’m not certain what it’s supposed to be.
Where is it suuposed to start and stop? Is the circumnavigation intentional?
Light
2023년 4월 3일
Star Strider
2023년 4월 3일
Just out of curiosity, I’ll look up the great circle calculations in the morning and see if I can make it work. I may be able to find the error.
Star Strider
2023년 4월 3일
Thanks. I pulled up the link, and I’ll go through it in a few minutes.
If you are doing everything in degrees, all the relevant trig functions should be in degrees. Since ‘psi’ is in degrees, the ‘X, Y, Z’ functions should also be in degrees. I changed those and re-ran the code, however while improved, it still ends up in Africa. I’ll look at that in a few minutes.
Light
2023년 4월 3일
As always, my pleasure!
Here is my analysis, and it appears to provide the correct flight path (that should provide a nice view of tha aurora borealis about midway through the flight for those on the left side of the airplane) —
% % % % % GREAT CIRCLE NAVIGATION
% % % great_circle.m
% % % 03-Mar-2023
% % % Star Strider
b1 = 32.78; % StartLat
l1 = -96.81; % StartLon
b2 = 50.07; % EndLat
l2 = 14.42; % EndLon
r = 6371; % Earth Radius (km)
[x1,y1,z1] = sp2ct(b1,l1);
[x2,y2,z2] = sp2ct(b2,l2);
psi = acosd(x1.*x2+y1.*y2+z1.*z2) % Cities Angular Distance (°)
Distance = 2*pi*r*(psi/360);
fprintf('Distance = %.2f km\n',Distance)
phiv = linspace(0,psi,250).'; % Angular Flight Path
[xv,yv,zv] = angdist(x1,y1,z1,x2,y2,z2,phiv,psi);
[lv,bv] = ct2sp(xv,yv,zv);
lv(lv<0) = lv(lv<0) + 360;
figure
load('topo.mat', 'topo', 'topomap1');
% load earth topography
whos topo topomap1;
% set coordinates of points
Longitude = 0:15:360; Latitude = -72:6:72;
% set background parameters
contour(0:359, -89:90, topo, [0,0], 'w');
axis equal
box on
set(gca, 'Color', [0 0 1], 'XLim', [0 360], 'YLim', [-90 90], 'XTick',...
[0: 180: 360], 'YTick', [-90: 90: 90]);
%set(gca, 'Color', [0 0 1], 'XLim', [0 360], 'YLim', [-90 90]);
grid; hold on
plot(lv, bv, 'y.', 'linewidth', 1);
hold off
function [lon,lat] = ct2sp(x,y,z)
lon = atan2d(y,x);
lat = asind(z);
end
function [x,y,z] = angdist(x1,y1,z1,x2,y2,z2,phi,psi)
% x90 = (x2-x1.*cosd(psi))./sind(psi);
% y90 = (y2-y1.*cosd(psi))./sind(psi);
% z90 = (z2-z1.*cosd(psi))./sind(psi);
% x = x1.*cosd(phi).*x90.*sind(phi);
% y = y1.*cosd(phi).*y90.*sind(phi);
% z = z1.*cosd(phi).*z90.*sind(phi);
x = (x1.*sind(psi-phi)+x2.*sind(phi))./sind(psi);
y = (y1.*sind(psi-phi)+y2.*sind(phi))./sind(psi);
z = (z1.*sind(psi-phi)+z2.*sind(phi))./sind(psi);
end
function [x,y,z] = sp2ct(lat,lon)
x = cosd(lon).*cosd(lat);
y = sind(lon).*cosd(lat);
z = sind(lat);
end
I created functions for the repetitive calculations, and also to simplify the code design. (This is simply a personal preference, and I do not intend it to be criticism of your code.) There are also MATLAB functions cart2sph and sph2cart that will do the conversions. I coded the versions of them here that correspond to the provided great circle documentation. I will let you sort this to determine where the errors may be in your code.
For my part, I now have a great circle .m file if I want to use it in the future!
.
Light
2023년 4월 3일
Star Strider
2023년 4월 3일
Thank you!
That was actually fun, although it took me a bit longer to code it than I had expected it to.
The functions here are required to be at the end of the calling code (the order in which the appear is not important, and they can be called anywhere within it and by each other if necessary, just as can other such functions), although anonymous functions (these are not such) must appear in the code before they are used.
The plotting section has to be the last part of the calling code, after all the variables it uses are calculated, and before the function section. Note that the functions are all vectorised, and no (explicit) loops are required.
.
Light
2023년 4월 4일
Star Strider
2023년 4월 4일
There are three functions: ‘ct2sp’, ‘angdist’, and ‘sp2ct’. I named the first and the last so as not to overshadow similar core MATLAB functions.
All the calculations are from the documentation you provided. The ‘angdist’ function is from 
, although I initially coded it as 
. When I read down that far, it was preferable because it was obviously more efficient. The earlier code I kept in, although commented-out so that it would not execute.
, although I initially coded it as . When I read down that far, it was preferable because it was obviously more efficient. The earlier code I kept in, although commented-out so that it would not execute. The φ value is actually a vector:
phiv = linspace(0,psi,250).'; % Angular Flight Path
and follows the observation:
‘If φ
, then you are in the first city. If φ=ψ, then you are in the second city.’
, then you are in the first city. If φ=ψ, then you are in the second city.’ so that extends from the departure airport (0) to the destination airport (ψ), where ψ was defined as:
psi = acosd(x1.*x2+y1.*y2+z1.*z2) % Cities Angular Distance (°)
All the calculations are in degrees.
.
Light
2023년 4월 4일
Star Strider
2023년 4월 4일
‘Okay,its mainly following the guide given by the calculations guide.’
Pretty much completely following it.
‘Could you just explain the the 'angdist' function and the 'function' subplot at the end?’
I thought I already explained ‘angdist’ as being 
applied to each ‘(x1,x2)’, ‘(y1,y2)’, and ‘(z1,z2)’ pair in turn. It is the same relation applied to each pair. It returns the Cartesian coordinates ‘xv’, ‘yv’, and ‘zv’ that are vectors because φ is the vector ‘phiv’. Those are then transformed into the longitude coordinates ‘lv’, and corresponding latitude coordinates ‘bv’, in the ‘ct2sp’ call immediately after that, both of them being vectors as well, because their arguments are vectors. The ‘lv’ vector is further processed by adding 360 to its negative elements. Those are then plotted on the map.
applied to each ‘(x1,x2)’, ‘(y1,y2)’, and ‘(z1,z2)’ pair in turn. It is the same relation applied to each pair. It returns the Cartesian coordinates ‘xv’, ‘yv’, and ‘zv’ that are vectors because φ is the vector ‘phiv’. Those are then transformed into the longitude coordinates ‘lv’, and corresponding latitude coordinates ‘bv’, in the ‘ct2sp’ call immediately after that, both of them being vectors as well, because their arguments are vectors. The ‘lv’ vector is further processed by adding 360 to its negative elements. Those are then plotted on the map. There is no function subplot. The code after the figure call I copied from your existing code because it loads the map and sets its parameters. (I cannot improve on that, and besides it is consistent with the desired result with respect to displaying the flight path.) The only difference between my code and yours is:
plot(lv, bv, 'y.', 'linewidth', 1);
The three functions at the end of my script are the functions I use in the code.
I just now noticed that you did not post the MATLAB version/release that you are using in either of your threads. My code should run without error in all recent releases. If you are having problems with it because of version/release compatibility issues, I might be able to help with that. The online Run feature here and my offline computers all run R2023a.
.
Light
2023년 4월 4일
Star Strider
2023년 4월 4일
As always, my pleasure!
For that, see the documentation section on Create Functions in Files. This is relatively new (becoming available in the last few years). You will definitely have the ability to use them in R2022b.
I learned some things from answering your Question, so thank you for posting it!
Light
2023년 4월 4일
Star Strider
2023년 4월 4일
As always, my pleasure!
Light
2023년 4월 12일
Star Strider
2023년 4월 12일
It would, although I would have to see it before I comment further.
Light
2023년 4월 12일
편집: Walter Roberson
2023년 4월 21일
Today turned out to be an unscheduled ‘errand day’ so I was away for several hours.
If I am correct, the objective here is to develop a set of differnetial equations and then estimate parameters using the supplied data. If that is the situation, it would be relatively straightforward to code that. See Parameter Estimation for a System of Differential Equations for an example.
Extracting ‘VAR DATA’ to a table and saving it to a .mat file is straightforward.
VAR_DATA = array2table([0 108
0.5 107
1.1 106
1.7 105
2.3 104
2.9 104
3.5 103
4.0 102
4.6 101
5.2 100
5.7 99
6.3 99
6.8 98
7.4 97
7.9 97
8.4 96
9.0 95
9.5 94
10.0 94
10.6 93
11.1 92
11.6 92
12.1 91
12.6 91
13.1 90
13.6 89
14.1 89
14.6 88
15.1 88
15.6 87
16.1 87
16.5 86
17.0 86
17.5 85
18.0 85
18.4 84
18.9 83
19.4 82], 'VariableNames',{'Distance (m)','Speed (kph)'})
save('VAR_DATA.mat','VAR_DATA')
L = size(VAR_DATA,1);
t = linspace(0, L-1, L);
figure
stem3(t, VAR_DATA{:,1}, VAR_DATA{:,2}, 'filled')
grid on
xlabel('t')
ylabel('Distance')
zlabel('Velocity')
axis('equal')
figure
plot(VAR_DATA{:,1}, VAR_DATA{:,2}, '.-')
grid on
xlabel('Distance')
ylabel('Velocity')
axis('equal')
I will address this tomorrow. In the interim, using this information and my prototype code, experiment with coding it. Ideally, it would be preferable if there was also a time vector. in addition to distance and velocity. It it possible to get that, or infer it from the data? For example, are the VAR data sampled at a constant sampling interval or frequency? (I don’t know, and I’ve never encountered this previously, so I’ve never explored it. I just enjoy watching the matches.)
However, it might be possible to use the Symbolic Math Toolbox to solve the differential equations in a closed form , with that solution being the objective function for a nonlinear parameter estimation problem. Try that as well.
.
Light
2023년 4월 13일
Star Strider
2023년 4월 13일
No worries. I’m still catching up from yesterday. I’ll take another look at it and make some suggestions.
Light
2023년 4월 19일
편집: Walter Roberson
2023년 4월 21일
I created and uploaded ‘VAR_DATA.mat’ to make this a bit easier.
I am still somewhat lost on this, and have not given it much thought.
syms k t x(t) y(t) Y T
m = sym(450); % g
g = sym(9.81); % m/s^2
v0 = sym(108); % km/h
v = sqrt(diff(x)^2+diff(y)^2);
Eq1 = m*diff(x,2) == -k*v*diff(x);
Eq2 = m*diff(y,2) == -m*g - k*v*diff(y);
[VF,Subs] = odeToVectorField(Eq1,Eq2)
deqfcn = matlabFunction(VF,'Vars',{T,Y,k})
load('VAR_DATA.mat')
VAR_DATA
VAR_MTX = table2array(VAR_DATA);
t = linspace(0, 3, size(VAR_MTX,1)).';
B0 = [1; pi/2; VAR_MTX(1,2)];
B = lsqcurvefit(@objfcn, B0, t, VAR_MTX)
xy = objfcn(B,t)
figure
plot3(VAR_MTX(:,1), VAR_MTX(:,2), t, '.')
hold on
plot3(xy(:,1), xy(:,2), t, '-', 'LineWidth',1.5)
grid
xlabel('Distance')
ylabel('Velocity')
zlabel('t')
legend('Data','Fit', 'Location','best')
function xy = objfcn(b,t)
k = b(1);
theta = b(2);
v0 = b(3);
deqfcn = @(T,Y,k)[Y(2);k.*sqrt(Y(2).^2+Y(4).^2).*Y(2).*(-1.0./4.5e+2)-9.81e+2./1.0e+2;Y(4);k.*sqrt(Y(2).^2+Y(4).^2).*Y(4).*(-1.0./4.5e+2)];
ic = [v0*sin(theta) 0 v0*cos(theta) 0];
[t,y] = ode45(@(t,y)deqfcn(t,y,k), t, ic);
xy = y(:,[3 1]);
end
This is how I usually solve these sorts of problems, however I genuinely do not have a good feeling about what the variables are here. We do not have a time vector (although the exact time vector would be helpful), and I am not certain how to model this. I used the Symbolic Math Toolbox to derive the code for the differential equaations simply because it avoids the algebraic errors that I usually make (and that can be frustrating to detect and correct). Experiment with this code to estimate the parameters and correctly determine what the data are that are being modeled.
Note that you can select the variables that ‘objfcn’ returns by defining what columns the ‘xy’ variable selects from the four columns of ‘y’ returned by the ode45 call. These values will be used by lsqcurvefit to fit the data. The ‘names’ of the variables are in the ‘Subs’ output of the odeToVectorField call.
The earlier problem was straightforward aand relatively easy to solve. This one is not.
.
Light
2023년 4월 20일
Star Strider
2023년 4월 20일
As always, my pleasure!
My problem is that it should be obvious that either x or y is distance and either 
or 
is velocity, I cannot figure out how they relate here. I’ve been away from problems like this for longer than I’d like to admit, concentrating instead on chemical kinetics, physiological modeling, and signal processing.
or is velocity, I cannot figure out how they relate here. I’ve been away from problems like this for longer than I’d like to admit, concentrating instead on chemical kinetics, physiological modeling, and signal processing. It’s not obvious to me how to set this up, given the provided information.
Light
2023년 4월 20일
Star Strider
2023년 4월 20일
With respect to 3), the derivatives and ‘x’ vectors would be returned from the ode45 call (in my code, columns 2, 3, and 4). I would then compute the integral as trapz(x, sqrt(1 + dy./dx)).
If I could figure out how the results of integrating the differential equations figured into fitting the data, this would likely be straightforward. I still have no idea how to match those.
Light
2023년 4월 20일
Star Strider
2023년 4월 20일
If I could figure out what ‘x’ and ‘y’ corresponded to in the data, then yes. My code is designed to fit the data and estimate the parameters, however that requires knowing how to fit the integrated differential equation to the data, and that continues to elude me. (I’ve tried various conbinations of x and y to the position data, and various combinations of 
and 
to the velocity data, as well as various other combinations, and none of them have worked.)
and to the velocity data, as well as various other combinations, and none of them have worked.)
Light
2023년 4월 20일
g = ...;
k = ...;
m = ...;
[T,Y] = ode45(@(t,y)fun(t,y,g,k,m),[0 20],[1 1 1 1])
function dy = fun(t,y,g,k,m)
normv = sqrt(y(3)^2+y(4)^2);
dy = zeros(4,1);
dy(1) = y(3);
dy(2) = y(4);
dy(3) = -k*normv*y(3)/m;
dy(4) = (-m*g-k*normv*y(4))/m;
end
Light
2023년 4월 20일
Light
2023년 4월 20일
Light
2023년 4월 20일
m = 0.45; % in kg
g = 9.81; % gravitational acceleration in m/s^2
v0 = 108; % in km/h
v0 = v0*1000/3600; % in m/s
t = 0: 0.1: 2;
theta = 30*pi/180; %initial guess for theta
k = 3; %initial guess for air resistance in N
ICs = [0 0 v0*cos(theta) v0*sin(theta)];
[T,Y] = ode45(@(t,y)fun(t,y,g,k,m),t,ICs);
plot(T,Y(:,1:2))
function dy = fun(t,y,g,k,m)
normv = sqrt(y(3)^2+y(4)^2);
dy = zeros(4,1);
dy(1) = y(3);
dy(2) = y(4);
dy(3) = -k*normv*y(3)/m;
dy(4) = (-m*g-k*normv*y(4))/m;
end
Light
2023년 4월 20일
could it be written using the same notation like what I had in my last comment with:
mx2 = @(t,x) [x(1); x(2); -k*mag_v*x(1)]
my2 = @(t,y) [y(1); y(2); -m*g -k*mag_v*y(2)]
No. The 4 equations must be solved together setting the solution vector to z = [x, y, dx/dt, dy/dt].
And because mag_v has to be updated according to the values of dx/dt and dy/dt, the dz/dt vector should be set in a function, not as a function handle. Just as I did it.
Light
2023년 4월 20일
Light
2023년 4월 20일
Does mag_v need to be changed to normv? Or it can remain as mag_v?
You can name a variable as you like it.
Also how can you see the values for X and Y? From the plot seen above
X(t) is saved in Y(:,1), Y(t) is saved in Y(:,2) as you can see in the plot command:
plot(T,Y(:,1:2))
To see the values, just add the lines
Y(:,1)
Y(:,2)
It seems it's time that you learn MATLAB fundamentals:
Light
2023년 4월 21일
Light
2023년 4월 21일
Additionally solve the ODE
ds/dt = sqrt(1+(y(4)/y(3))^2)
with inital condition
s(0) = 0.
Thus your y-vector now becomes 5x1 instead of 4x1 with s being the solution component y(5).
Light
2023년 4월 21일
Star Strider
2023년 4월 21일
With respect to using trapz, use cumtrapz instead to get a vector result rather than a scalar result.
Light
2023년 4월 21일
Star Strider
2023년 4월 21일
Yes. The vector lengths have to be the same, so rather than the colon operator, linspace would be best, unless the ‘x1’ vector was already defined as something else. If it is, use that vector instead.
Star Strider
2023년 4월 21일
The ‘f’ value is a scalar, so cumtrapz takes that as the dimension argument and returns the error.
I am sort of lost at this point, however if ‘f’ is supposed to be a vector and if the ‘y’ values you want to use are the integrated results of the ode45 call, consider defining ‘f’ as:
f = sqrt(1 + (y(:,4) ./ y(:,3)).^2);
It will be a column vector as well, since its arguments are column vectors.
The interp1 result needs to be saved to a variable. (It will return the value of ‘Y(:,2)’ corresponding to the interpolated value of ‘Y(:,1)’ that equals 0.7.)
I would also leave the values of ‘x’ and ‘y’ as column vectors for consistency, rather than transposing them. Any row vectors created (for example from linspace) need to be transposed to column vectors as well, again for consistency.
.
Light
2023년 4월 21일
Am I really that unclear in what I explain to you ?
m = 0.45; % in kg
g = 9.81; % gravitational acceleration in m/s^2
v0 = 108; % in km/h
v0 = v0*1000/3600; % in m/s
t = 0: 0.01: 2;
theta = 30*pi/180; %initial guess for theta
k = 3; %initial guess for air resistance in N
ICs = [0 0 v0*cos(theta) v0*sin(theta) 0]; % <- Add a 0 at position 5 for s(0)=0
[T,Y] = ode45(@(t,y)fun(t,y,g,k,m),t,ICs);
% Plot s
plot(T,Y(:,5)) % <- Plot s
function dy = fun(t,y,g,k,m)
normv = sqrt(y(3)^2+y(4)^2);
dy = zeros(4,1);
dy(1) = y(3);
dy(2) = y(4);
dy(3) = -k*normv*y(3)/m;
dy(4) = (-m*g-k*normv*y(4))/m;
dy(5) = sqrt(y(3)^2+y(4)^2); %<- Add the (correct) differential equation for s
end
Withiout running the complete code, so creating a random matrix for ‘Y’ —
Y = randn(10,4)
x = Y(:,1);
y = Y(:,2);
f = sqrt(1 + (Y(:,4) ./ Y(:,3)).^2)
s = cumtrapz(x, f)
So perhaps something like this?
.
Light
2023년 4월 21일
Torsten
2023년 4월 21일
As noted in your other post, the equation for s is
ds/dt = sqrt((dx/dt)^2 + (dy/dt)^2)
not
ds/dt = sqrt(1+ ((dy/dt)/(dx/dt))^2).
It already follows from the unit of s (meter) that the second version cannot be correct.
If it were, s had unit "second".
Light
2023년 4월 21일
If s is the distance travelled until t=x, the correct formula is
s(x) = integral_{t=0}^{t=x} sqrt((dx/dt)^2 + (dy/dt)^2) dt
not
s(x) = integral_{t=0}^{t=x} sqrt(1 + ((dy/dt)/(dx/dt))^2) dt
Light
2023년 4월 21일
Light
2023년 4월 21일
추가 답변 (0개)
카테고리
도움말 센터 및 File Exchange에서 Lengths and Angles에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
