Vectorize my code in azimuth calculation.

6 次查看（过去 30 天）
Mani Ahmadian
Mani Ahmadian 2014-9-23
评论： Mani Ahmadian 2014-9-25
Hi I have lots of points(in Latitude & Longitude) and I want to calculate the azimuth of each pairs of them. I use code as attached file 'AzimuthCalc.m'. when I run my function, it's so time consuming and it's elapsed time is 109.371341 seconds. I'm sure the problem caused by my non-vectorize code.
So, would you please check my function and help me to improve it?
Thanks a lot Mani
  2 个评论
Adam
Adam 2014-9-23
Run:
profile on
your code
profile off
profile viewer
and find out where your code is slow. No point in guessing when you can use facts!
Mani Ahmadian
Mani Ahmadian 2014-9-23
Dear Adam,
Thanks so much. I guess the problem is here(nested for):
for ii=1:length(lat)
for jj=ii+1:length(lat)
AzimuthPoints(ii,jj)=azimuth('rh',lat(ii,1),long(ii,1),lat(jj,1),long(jj,1));
end
end
If I can change my simple code to vectorize one, problem will solve. So, I'm trying to find a way to do that, so plz help me.
Thanks
Mnai

请先登录，再进行评论。

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

采纳的回答

Guillaume
Guillaume 2014-9-23
I don't have the mapping toolbox, so I'm just relying on the documentation to write my answer.
For a start, you can replace the inner loop with:
AzimuthPoints(ii,ii+1:end)=azimuth('rh',lat(ii),long(ii),lat(ii+1:end),long(ii+1:end));
I believe the following may work to replace both loops (and the AzimuthPoints declaration)
[lat1, lat2] = meshgrid(lat);
[long1, long2] = meshgrid(long);
AzimuthPoints = azimuth('rh', lat1, long1, lat2, long2);
The downside of this method is that it's going to calculate the whole AzimuthPoints matrix, not just the upper triangle, so it may actually be slower.
You could do:
lat1 = triu(lat1, 1);
%same with lat2, long1 and long2
to replace the lower triangle and main diagonal with zero, on the assumption that calculating
azimuth('rh', 0, 0, 0, 0)
is faster.
  5 个评论
显示 3更早的评论
Guillaume
Guillaume 2014-9-24
TrueAzimuth should be 0s and 1s but of type logical. If you used my code to generate it, it shoul already be logical, but if not:
TrueAzimuth = logical(TrueAzimuth);

请先登录，再进行评论。

更多回答（1 个）

Roger Stafford
Roger Stafford 2014-9-24
You could try this. It calls on 'azimuth' in the same sequence as in your for-loop code, but it is vectorized. Whether it is faster depends on whether the 'tril' and 'repmat' calls take less time than the overhead in the for-loops and individual 'azimuth' calls in your code.
n = length(lat);
p = tril(repmat(1:n,n,1),-1);
p = p(p>0);
q = tril(repmat((1:n)',1,n),-1);
q = q(q>0);
AzimuthPoints = zeros(n);
AzimuthPoints(p+n*(q-1)) = azimuth('rh',lat(p),long(p),lat(q),long(q));

类别

在 Help Center 和 File Exchange 中查找有关 Time Series Objects 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by