Plotting Right Ascension & Declination on Mercator Map

I'm developing a program to do many different satellite orbit calculations. One of them involves plotting the ground tracks of a given satellite on a mercator (lat & long) map. I have the Right Ascension and Declination angles in a vector called "rad" in the following code. With Right Ascension in the first column, and declination in the second column. The plot is coming out okay, except that it is connecting the groundtrack line when the right ascension (width) resets from 360 to 0. I can't quite figure out how to plot the data without connecting the line when the right ascension resets to 0 degrees from 360.

rad(:,1) ranges from 0 to 360 (right ascension) rad(:,2) ranges from -90 to 90 (declination)

                earth = imread('ImageOfEARTH');

                hold on

                F = earth(end:-1:1,:,:);

                image([0 360],[-90 90],F)

                plot(rad(:,1),rad(:,2),'red');

                axis([0,360,-90,90])
                daspect([1 1 1]);

This code outputs the following plot:

See how that line is connected through the middle? I need to eliminate that.

Thanks!

-Brandon

Answer 1

[Star Strider]

I don’t have ‘ImageOfEARTH’ so I can’t run your code. However, one way to deal with the ‘wrap-around’ problem is to add NaN at the end of each of your vectors. Those will ‘break’ the continuity of your plot, since NaN does not plot, and your line plot should no longer wrap.

[Brandon Johnson]

Try this for the image:

https://upload.wikimedia.org/wikipedia/commons/thumb/e/ea/Equirectangular-projection.jpg/800px-Equirectangular-projection.jpg

And as for rad:

rad =

212.5071   20.5257
211.5328   18.4516
210.5980   16.4317
209.6972   14.4616
208.8256   12.5371
207.9791   10.6545
207.1542    8.8102
206.3478    7.0011
205.5572    5.2239
204.7799    3.4759
204.0137    1.7542
203.2565    0.0565
202.5066   -1.6197
201.7621   -3.2767
201.0214   -4.9166
200.2829   -6.5416
199.5450   -8.1536
198.8063   -9.7545
198.0652  -11.3462
197.3201  -12.9304
196.5696  -14.5090
195.8121  -16.0837
195.0457  -17.6560
194.2688  -19.2278
193.4795  -20.8007
192.6757  -22.3762
191.8552  -23.9562
191.0155  -25.5421
190.1540  -27.1356
189.2675  -28.7385
188.3527  -30.3524
187.4058  -31.9789
186.4222  -33.6199
185.3970  -35.2770
184.3242  -36.9520
183.1970  -38.6464
182.0071  -40.3621
180.7450  -42.1007
179.3992  -43.8636
177.9558  -45.6523
176.3977  -47.4678
174.7042  -49.3109
172.8493  -51.1816
170.8008  -53.0793
168.5175  -55.0020
165.9472  -56.9457
163.0223  -58.9040
159.6558  -60.8659
155.7344  -62.8146
151.1130  -64.7241
145.6098  -66.5548
139.0109  -68.2476
131.0982  -69.7173
121.7225  -70.8469
110.9397  -71.4917
 99.1538  -71.5007
 87.1171  -70.7557
 75.6755  -69.2055
 65.4323  -66.8674
 56.6120  -63.7965
 49.1508  -60.0531
 42.8496  -55.6836
 37.4804  -50.7163
 32.8354  -45.1646
 28.7400  -39.0338
 25.0517  -32.3289
 21.6533  -25.0633
 18.4444  -17.2677
 15.3348   -8.9995
 12.2363   -0.3480
  9.0552    8.5623
  5.6826   17.5772
  1.9803   26.5223
357.7613   35.2158
352.7573   43.4785
346.5696   51.1325
338.5985   57.9875
327.9954   63.8037
313.8302   68.2489
295.9686   70.9159
276.4840   71.5519
258.8962   70.3807
245.1198   67.9986
234.9275   64.9726
227.3835   61.6688
221.6610   58.2895
217.1848   54.9394
213.5778   51.6709
210.5925   48.5085
208.0641   45.4615
205.8793   42.5312
203.9589   39.7145
202.2458   37.0059
200.6979   34.3992
199.2836   31.8875
197.9787   29.4640
196.7646   27.1221
195.6261   24.8556
194.5515   22.6587
193.5310   20.5257

[Brandon Johnson]

And for the NaN? Should I place it entirely in the third column?

[Star Strider]

Assuming all your ‘rad’ data are similar in their construction for all your satellites (that is, the wrap-around occurs in the interior and not at the end of the matrix), this works:

[~,nan_ins] = min(rad(:,1));                        % Index Of Mimimum ‘rad(:,1)’
rad_nan = rad(1:nan_ins,:);                         % Create First Section Of ‘rad_nan’
rad_nan(nan_ins+1,:) = [NaN NaN];                   % Insert ‘NaN’ Row
rad_nan = [rad_nan; rad(nan_ins+1:end,:)];          % Append Remaining Rows Of ‘rad’ To ‘rad_nan’

although it looks a bit kludgy.

Remember to plot:

plot(rad_nan(:,1),rad_nan(:,2),'red');

instead of ‘rad’.

[Brandon Johnson]

Beautiful, now I get this plot. But the problem is that now if I try to plot multiple periods of the satellite (multiple times around) I will have many discontinuities of this type. Is there a way to modify that code to handle multiple discontinuities?

[Star Strider]

First, thanks for the earth map link. The image is now in one of my MATLAB directories.

Second, it depends on how your data are stored and the lengths of the individual (Rx2) matrices for each satellite revolution. If there are N of them and if they are all the same lengths (that is 100 rows each so R=100), then you can use the reshape function to convert the long (R*Nx2) matrix into N (Rx2) matrices in the same matrix (so it would be (Rx2*N)), and then use the same code to insert a (2*N) row of NaN values in it.

Otherwise, you will have to test for the wrap-around (minimum) for each two-column section in your (Rx2) matrix. If the minima vary but are all <2, then use the find function to test for rad(:,1)<2, or whatever test is appropriate. You would likely have to use a loop to use my code with each section of the matrix. So if possible, make them all R rows long, use reshape, and avoid the complexity.

Third, if all the (Rx2) matrix segments ‘break’ at the same row, you can simplify my code a bit by using the circshift function to rotate your RxN*2 reshaped matrix so that all the vectors in your RxN*2 reshaped matrix begin at the same row. Then you simply add a 2*N row of NaN values to the end of the matrix and be done with it. That’s the easiest way. (I thought about doing it with your matrix here.) You would then use reshape again to create one long (N*(R+N)x2) matrix, now including the added NaN rows.

Also, the satellite track has to begin somewhere and the calculation end somewhere else, so some discontinuity is inevitable (unless it ends exactly over its launch site).

Answer 2

[Chad Greene]

As Star Strider mentioned, you can add a NaN where it wraps around. A faster and possibly less visually appealing way is to plot red dots instead of a red line:

plot(rad(:,1),rad(:,2),'r.');

A note unrelated to your question: Your map isn't actually a Mercator map. It's simply unprojected geographic coordinates where you've given all degrees longitude and latitude equal size. Mercator projections make Greenland appear larger than Africa.

[Brandon Johnson]

Yeah, I know the map isn't right. I just used it for sake of demonstration.

So you're saying I should write a routine to find the point where it wraps around (or multiple points) and then add the NaN to the third element? Or simply replace the value with NaN?

[Brandon Johnson]

I've replaced a row with NaNs and I'm getting the same result.