McGill University
School of Computer Science

Computer Science COMP 199 (Winter term)
Excursions in Computing Science

%MATLABpak_ii: worldPopMyAnc.m, Galton.m, drawStars.m, selectStarsUpto52ly, drawConst3D.m

% worldPopMyAnc.m		THM		080729
% Plot of world population, linear and semilog.
% Data from en.wikipedia.org/wiki/World_population
%X=[-70000,-10000,-9000,-8000,-7000,-6000,-5000,-4000,-3000,-2000,-1000,-500,1,1000,1750,1800,1850,1900,1950,1955,1960,1965,1970,1975,1980,1985,1990,1995,2000,2005];
X=[-10000,-9000,-8000,-7000,-6000,-5000,-4000,-3000,-2000,-1000,-500,1,1000,1750,1800,1850,1900,1950,1955,1960,1965,1970,1975,1980,1985,1990,1995,2000,2005];
%Y=[2,1000,3000,5000,7000,10000,15000,20000,25000,35000,50000,100000,200000,310000,791000,978000,1262000,1650000,2518629,2755823,2981659,3334874,3692492,4068109,4434682,4830979,5263593,5674380,6070581,6453628];
Y=[1000,3000,5000,7000,10000,15000,20000,25000,35000,50000,100000,200000,310000,791000,978000,1262000,1650000,2518629,2755823,2981659,3334874,3692492,4068109,4434682,4830979,5263593,5674380,6070581,6453628];
hold off
plot(X,Y/1000000)
hold on
Xanc=[2000:-25:1175];
Yanc=[2.^((2000.-Xanc)./25)];
plot(Xanc,Yanc/1e9,'g')
xlabel('Year (BCE to present)')
ylabel('World Population / My Ancestors (billions)')
hold off
% OR
semilogy(X,Y*1000)
hold on
semilogy(Xanc,Yanc,'g')
xlabel('Year (BCE to present)')
ylabel('World Population / My Ancestors')
hold off


% function distrib = Galton(n,h)		THM			090916
% For a Galton board of h levels simulate resulting distribution for n balls:
% use h tests rand(1)>0.5 and accumulate in histogram cell whose index os the
% sum of the number of right-bounces:
% e.g.	0	1	2	3	4
%	llll	lllr	llrr	lrrr	rrrr
%		llrl	lrlr	rlrr
%		lrll	rllr	rrlr
%		rlll	lrrl	rrrl
%			rlrl
%			rrll
function distrib = Galton(n,h)
distrib = zeros(1,h+1);
for j = 1:n
 right = 0;
 for k = 1:h
  right = right + (rand(1)>0.5);	% need the parentheses
 end % for k
 distrib(right+1) = distrib(right+1)+1;
end % for j


% function drawStars(stars)
% THM 070502	in file drawStars.m	Copied 080731 from lecture1/matlab
% stars is n*19 array of doubles (e.g., created by selectStars)
function drawStars(stars)
XYZmax = -Inf; magmax = -Inf; magmin = Inf; parxmax = -Inf; parxmin = Inf;
siz = size(stars);
for k = 1:siz(1)
   RA = pi/180*(360/24*(stars(k,9) + 1/60*(stars(k,10) + 1/60*stars(k,11))));
									% radian
   sD = (stars(k,12)=='+')*1 + (stars(k,12)=='-')*(-1);			% sign
   DEC = pi/180*sD*(stars(k,13) + 1/60*(stars(k,14) + 1/60*stars(k,15)));
									% radian
   ly = 3.26/stars(k,19);						% dist
   %ly = 1;						% show on unit sphere
   xm = ly/105;							% exameters
   X(k) = xm*cos(DEC)*cos(RA);
   Y(k) = xm*cos(DEC)*sin(RA);
   Z(k) = xm*sin(DEC);
   if XYZmax < X(k), XYZmax = X(k); end			% max for plot axes
   if XYZmax < Y(k), XYZmax = Y(k); end
   if XYZmax < Z(k), XYZmax = Z(k); end
   if magmax < stars(k,16), magmax = stars(k,16); end
   if magmin > stars(k,16), magmin = stars(k,16); end
   if parxmax < stars(k,19), parxmax = stars(k,19); end
   if parxmin > stars(k,19), parxmin = stars(k,19); end
end
XYZmax, magmax, magmin, parxmax, parxmin
U = zeros(size(X));
V = zeros(size(Y));
W = zeros(size(Z));
axis([-XYZmax XYZmax -XYZmax XYZmax -XYZmax XYZmax]),grid
xlabel('x'), ylabel('y'), zlabel('z')
quiver3(X,Y,Z,U,V,W,0,'r*')
hold on
quiver3(U,V,W,U,V,W,0,'ro')		% include Sol
hold off

% function starsUpto52ly = selectStarsUpto52ly()
% THM 070502	in file selectStarsUpto52ly.m
%  Adapted 080731 from lecture1/matlab/selectStars.m  NB 52.3ly is 1/2 Exameter
function starsUpto52ly = selectStarsUpto52ly()
yaleBSC = textread('/Users/tim/linux/pedagogy/astron/yale_bsc/yale.txt','%s','delimiter','\n','whitespace','');
for k = 1:size(yaleBSC)
   StarNo(k) = str2double(yaleBSC{k}(1:4));	% number in the Yale BSC
   RA19h(k) = str2double(yaleBSC{k}(53:54));	% right ascension, 1900: h, m, s
   RA19m(k) = str2double(yaleBSC{k}(56:57));
   RA19s(k) = str2double(yaleBSC{k}(59:62));
   DEC19sgn(k) = yaleBSC{k}(64);		% declination, 1900: +-, d, m, s
   DEC19d(k) = str2double(yaleBSC{k}(65:66));
   DEC19m(k) = str2double(yaleBSC{k}(68:69));
   DEC19s(k) = str2double(yaleBSC{k}(71:72));
   RA20h(k) = str2double(yaleBSC{k}(74:75));	% right ascension, 2000: h, m, s
   RA20m(k) = str2double(yaleBSC{k}(77:78));
   RA20s(k) = str2double(yaleBSC{k}(80:83));
   DEC20sgn(k) = yaleBSC{k}(85);		% declination, 1900: +-, d, m, s
   DEC20d(k) = str2double(yaleBSC{k}(86:87));
   DEC20m(k) = str2double(yaleBSC{k}(89:90));
   DEC20s(k) = str2double(yaleBSC{k}(92:93));
   vMag(k) = str2double(yaleBSC{k}(109:113));	% visual magnitude("brightness")
   annMoRA(k) = str2double(yaleBSC{k}(156:161));% annual motion of rt ascension
   annMoDEC(k) = str2double(yaleBSC{k}(163:168));% annual motion of declination
   parx(k) = str2double(yaleBSC{k}(171:174));	% parallax (seconds)
   rankB(k) = k;				% rank for sorting brightness
   rankC(k) = k;
end
nearank = sortrows([parx',rankC'],-1);		% descending by parallax
for k = 1:size(yaleBSC)
   ck = nearank(k,2);				% select in order of closeness
   if 3.26/parx(ck) > 52.3 break, end
   starsUpto52ly(k,1) = StarNo(ck);
   starsUpto52ly(k,2) = RA19h(ck);
   starsUpto52ly(k,3) = RA19m(ck);
   starsUpto52ly(k,4) = RA19s(ck);
   starsUpto52ly(k,5) = DEC19sgn(ck);
   starsUpto52ly(k,6) = DEC19d(ck);
   starsUpto52ly(k,7) = DEC19m(ck);
   starsUpto52ly(k,8) = DEC19s(ck);
   starsUpto52ly(k,9) = RA20h(ck);
   starsUpto52ly(k,10) = RA20m(ck);
   starsUpto52ly(k,11) = RA20s(ck);
   starsUpto52ly(k,12) = DEC20sgn(ck);
   starsUpto52ly(k,13) = DEC20d(ck);
   starsUpto52ly(k,14) = DEC20m(ck);
   starsUpto52ly(k,15) = DEC20s(ck);
   starsUpto52ly(k,16) = vMag(ck);
   starsUpto52ly(k,17) = annMoRA(ck);
   starsUpto52ly(k,18) = annMoDEC(ck);
   starsUpto52ly(k,19) = parx(ck);
end

% 080813
% size(startsUpto52ly) 178    19	NB I mis-typed it
% 
%   Columns 1 through 6
%   Columns 7 through 12 
%   Columns 13 through 18 
%   Column 19 
%  format long		NB last 11 places all 0
%  startsUpto52ly(:,19)
% StarNo RA19h RA19m RA19s DEC19sgn DEC19d
% DEC19m DEC19s RA20h RA20m RA20s DEC20sgn
% DEC20d DEC20m DEC20s vMag annMoRA annMoDEC
% parx
% 
% startsUpto52ly	NB 43 is '+', 45 is '-'
%        --------1900-------- ---------2000--------
%                 vv                    vv NB 43 is '+', 45 is '-'
%        h  m  s  pm  d  m  s  h  m  s  pm  d  m  s vMag aMra aMdc  parx
%  5459 14 32 483 45 60 25 22 14 39 362 45 60 50 07 -000 -036  007 0.750
%  5460 14 32 483 45 60 25 22 14 39 362 45 60 50 07  013 -036  007 0.750
%  2491 06 40 446 45 16 34 44 06 45 089 45 16 42 58 -015 -005 -012 0.378
%  1084 03 28 131 45 09 47 48 03 32 558 45 09 27 30  037 -010  000 0.304
%  8085 21 02 248 43 38 15 27 21 06 537 43 38 44 57  052  041  032 0.294
%  8086 21 02 263 43 38 15 14 21 06 552 43 38 44 30  060  041  031 0.294
%  2943 07 34 040 43 05 28 53 07 39 181 43 05 13 30  004 -007 -010 0.292
%  8387 21 55 426 45 57 11 48 22 03 213 45 56 47 10  047  039 -026 0.290
%  0509 01 39 253 45 16 27 51 01 44 040 45 15 56 15  035 -017  009 0.287
%  1325 04 10 402 45 07 48 30 04 15 163 45 07 39 10  044 -022 -034 0.209
%  7557 19 45 542 43 08 36 15 19 50 469 43 08 52 06  008  005  004 0.202
%  6752 18  0 241 43 02 31 22 18 05 272 43 02 29 58  040  003 -011 0.201
%  6401 17 09 117 45 26 27 21 17 15 209 45 26 36 10  053 -005 -011 0.188
%  0487 01 36 005 45 56 42 08 01 39 477 45 26 36 05  053 -005 -011 0.188
%  6402 17 09 118 45 26 27 17 17 15 207 45 36 06 04  053  004 -016 0.179
%  7703 20 04 377 45 36 21 09 20 11 118 45 69 39 40  047  006 -017 0.177
%  7462 19 32 332 43 69 29 28 19 32 215 43 57 48 58  034  011 -005 0.176
%  0219  0 43 030 43 57 17 06  0 49 060 43 66 10 56  036  012 -011 0.176
%  7665 19 58 549 45 66 26 12 20 08 432 45 21 24 56  057  010 -017 0.172
%  5568 14 51 374 45 20 57 49 14 57 279 45 43 04 11  043  030  007 0.162
%  1008 03 15 560 45 43 27 08 03 19 556 45 77 15 16  028  022  003 0.159
%  0098  0 20 300 45 77 49 03  0 25 453 45 19 06 04  046  001 -001 0.156
%  5544 14 46 464 43 19 30 57 14 51 232 43 56 11 53  059  003 -000 0.152
%  0486 01 35 599 45 56 42 15 01 39 474 45 56 11 41  058  003  000 0.152
%  8728 22 52 076 45 30 09 08 22 57 390 45 29 37 20  012  003 -002 0.149
%  8832 23 08 278 43 56 36 58 23 13 169 43 57 10 07  056  021  003 0.146
%  0222  0 43 082 43 04 46  0  0 48 229 43 05 16 50  057  008 -011 0.143
%  0077  0 14 517 45 65 27 45  0 20 042 45 64 52 30  042  017  012 0.141
%  6426 17 12 087 45 34 52 40 17 18 571 45 34 59 23  059  012 -002 0.141
%  1543 04 44 246 43 06 47 12 04 49 503 43 06 57 41  032  005  000 0.137
%  4374 11 12 509 43 32 05 31 11 18 109 43 31 31 45  049 -004 -006 0.137
%  4375 11 12 509 43 32 05 31 11 18 109 43 31 31 45  044 -004 -006 0.137
%  6623 17 42 326 43 27 46 45 17 46 275 43 27 43 15  034 -003 -007 0.133
%  7001 18 33 331 43 38 41 26 18 36 562 43 38 47 01  000  002  003 0.133
%  6416 17 11 277 45 46 31 55 17 19 031 45 46 38 02  055  010  002 0.131
%  0321 01 01 368 43 54 25 47 01 08 163 43 54 55 14  052  034 -016 0.130
%  8721 22 50 503 45 32 05 41 22 56 239 45 31 33 57  065  003 -002 0.130
%  0753 02 30 357 43 06 24 35 02 36 048 43 06 53 13  058  018  015 0.129
%  1982 05 40 164 45 22 27 17 05 44 264 45 22 25 19  062 -003 -004 0.128
%  1983 05 40 176 45 22 28 51 05 44 278 45 22 26 54  036 -003 -004 0.128
%  6927 18 22 515 43 72 41 22 18 21 032 43 72 43 58  036  005 -004 0.128
%  4983 13 07 124 43 28 23 06 13 11 523 43 27 52 41  043 -008  009 0.124
%  7722 20 09 031 45 27 19 53 20 15 173 45 27 01 58  057  012 -002 0.124
%  0857 02 47 429 45 13 10 30 02 52 320 45 12 46 10  060  004 -002 0.121
%  0493 01 37 039 43 19 46 57 01 42 297 43 20 16 07  052 -003 -007 0.119
%  4496 11 35 471 43 34 46  0 11 41 029 43 34 12 05  053 -000 -004 0.119
%  8181 21 18 108 45 65 49 07 21 26 267 45 65 21 59  042  001  008 0.118
%  2261 06 13 132 45 74 43 07 06 10 146 45 74 45 11  051  001 -002 0.117
%  4785 12 28 596 43 41 54 03 12 33 445 43 41 21 27  043 -007  003 0.117
%  5019 13 13 103 45 17 45 18 13 18 242 45 18 18 41  047 -011 -011 0.117
%  1136 03 38 274 45 10 06 06 03 43 148 45 09 45 48  035 -001  007 0.113
%  1614 04 55 511 45 05 52 16 05  0 490 45 05 45 12  062  006 -011 0.110
%  3815 09 29 397 43 36 15 45 09 35 395 43 35 48 36  054 -007 -002 0.109
%  0996 03 14 069 43 03  0 13 03 19 216 43 03 22 13  048  003  001 0.108
%  5235 13 49 553 43 18 53 56 13 54 410 43 18 23 52  027 -001 -004 0.108
%  0511 01 40 305 43 63 21 32 01 47 448 43 63 51 09  056  006 -002 0.107
%  4550 11 47 130 43 38 26 10 11 52 587 43 37 43 08  065  040 -058 0.107
%  4458 11 29 378 45 32 18 05 11 34 294 45 32 49 53  060 -007  008 0.106
%  2047 05 48 276 43 20 15 28 05 54 229 43 20 16 34  044 -002 -001 0.104
%  4540 11 45 291 43 02 19 42 11 50 416 43 01 45 53  036  007 -003 0.104
%  0939 03 02 028 45 72 17 35 03 02 155 45 71 54 09  055  000  000 0.103
%  6212 16 37 309 43 31 47 01 16 41 171 43 31 36 10  028 -005  004 0.102
%  0166  0 34 095 43 20 42 40  0 39 217 43 21 15 02  059 -005 -004 0.099
%  4523 11 41 449 45 39 57 22 11 46 310 45 40 30 01  049 -015  004 0.099
%  4825 12 36 355 45  0 54 03 12 41 395 45 01 26 58  037 -006  000 0.099
%  4826 12 36 355 45  0 54 03 12 41 395 45 01 26 58  037 -006  000 0.099
%  7936 20 40 105 45 25 37 49 20 46 056 45 25 16 16  041 -000 -002 0.098
%  1006 03 15 358 45 62 57 28 03 17 461 45 62 34 32  055  013  007 0.097
%  5340 14 11 060 43 19 42 11 14 15 396 43 19 10 57 -000 -011 -020 0.097
%  9038 23 47 317 43 74 59 13 23 52 250 43 75 32 41  064  003  001 0.097
%  2990 07 39 118 43 28 16 04 07 45 189 43 28 01 34  011 -006 -001 0.094
%  5868 15 41 353 43 07 39 59 15 46 265 43 07 21 11  044 -002 -001 0.094
%  6098 16 17 425 45 69 51 32 16 28 281 45 70 05 04  049  002  001 0.094
%  0660 02 10 568 43 33 46  0 02 17 032 43 34 13 28  049  012 -002 0.093
%  6171 16 31 062 45 02 06 40 16 36 213 45 02 19 29  057  005 -003 0.093
%  0937 03 01 508 43 49 13 53 03 09 040 43 49 36 48  040  013 -001 0.092
%  6806 18 06 187 43 38 27 05 18 09 374 43 38 27 27  064 -003 -005 0.092
%  1010 03 16 023 45 62 53 17 03 18 128 45 62 30 23  052  013  007 0.091
%  0637 02 06 237 45 51 18 52 02 10 256 45 50 49 28  061  021  007 0.090
%  0683 02 14 296 45 26 25 08 02 18 584 45 25 56 44  063 -002  005 0.089
%  5553 14 48 514 43 19 33 19 14 53 236 43 19 09 10  060 -005  002 0.087
%  8322 21 41 313 45 16 34 52 21 47 023 45 16 07 38  029  003 -003 0.087
%  3384 08 28 572 45 31 10 52 08 32 514 45 31 30 03  064 -011  008 0.086
%  4102 10 22 247 45 73 31 22 10 24 237 45 74 01 54  040 -000 -000 0.086
%  5618 15  0 294 43 48 02 36 15 03 473 43 47 39 16  048 -004  000 0.086
%  1674 05 03 477 45 57 36 33 05 05 306 45 57 28 22  047 -000  001 0.085
%  3538 08 49 223 45 05 03 21 08 54 178 45 05 26 04  060 -004  000 0.085
%  5933 15 51 500 43 15 59 16 15 56 271 43 15 39 42  039  003 -013 0.084
%  5897 15 46 196 45 63 07 19 15 55 084 45 63 25 50  029 -002 -004 0.083
%  8501 22 11 423 45 54 06 31 22 18 154 45 53 37 40  054  004 -007 0.083
%  9088 23 56 567 43 26 33 10  0 02 101 43 27 04 56  057  008 -010 0.083
%  0799 02 37 219 43 48 48 20 02 44 119 43 49 13 43  041  003 -001 0.082
%  3259 08 13 391 45 12 17 36 08 18 238 45 12 37 55  060  003 -010 0.082
%  4534 11 43 575 43 15 07 52 11 49 035 43 14 34 19  021 -005 -001 0.082
%  8430 22 02 212 43 24 51 24 22 07 006 43 25 20 42  038  003  000 0.082
%  0483 01 35 416 43 42 06 42 01 41 471 43 42 36 49  050  008 -001 0.081
%  1925 05 33 134 43 53 26 26 05 41 203 43 53 28 52  062  000 -005 0.081
%  1708 05 09 180 43 45 53 47 05 16 413 43 45 59 53  001  001 -004 0.080
%  3862 09 37 436 45 23 28  0 09 42 144 45 23 54 56  049 -004  003 0.080
%  4587 11 55 365 45 09 52 34 12  0 444 45 10 26 46  056  001 -005 0.080
%  1532 04 43 068 45 17 07 03 04 47 362 45 16 56 04  055  001  002 0.079
%  5011 13 11 487 43 09 56 48 13 16 464 43 09 25 27  052 -003  002 0.079
%  6518 17 25 184 43 67 23 26 17 25   0 43 67 18 23  064 -005  000 0.078
%  0695 02 17 580 45 24 16 14 02 22 325 45 23 48 59  052  002 -001 0.077
%  4112 10 24 138 43 56 29 36 10 30 375 43 55 58 50  048 -002 -000 0.077
%  6585 17 36 122 45 51 46 50 17 44 086 45 51 50 03  052 -000 -002 0.077
%  3079 07 48 315 45 34 27 14 07 52 156 45 34 42 19  050 -002  002 0.077
%  7602 19 50 240 43 06 09 25 19 55 187 43 06 24 24  037  000 -005 0.076
%  0159  0 32 125 45 25 19 03  0 37 206 45 24 46 02  056  014 -000 0.076
%  7957 20 43 153 43 61 27 01 20 45 173 43 61 50 20  034  001  008 0.075
%  0963 03 07 493 45 29 22 52 03 12 042 45 28 59 14  039  003  006 0.075
%  3569 08 52 218 43 48 26 04 08 59 124 43 48 02 30  031 -004 -002 0.075
%  5864 15 40 593 45 37 35 58 15 47 288 45 37 54 59  060 -004 -002 0.075
%  6458 17 16 550 43 32 35 47 17 20 395 43 32 28 04  054  001 -010 0.075
%  7578 19 48 184 45 24 11 23 19 54 176 45 23 56 28  062 -001 -004 0.075
%  8969 23 34 483 43 05 05 03 23 39 570 43 05 37 35  041  004 -004 0.074
%  0553 01 49 068 43 20 19 09 01 54 383 43 20 48 29  026  001 -001 0.074
%  3522 08 46 385 43 28 42 46 08 52 358 43 28 19 51  060 -005 -002 0.074
%  8729 22 52 331 43 20 13 57 22 57 279 43 20 46 08  055  002  001 0.074
%  3391 08 30 191 43 65 22  0 08 39 117 43 65 01 15  056 -000  001 0.073
%  3881 09 42 085 43 46 29 13 09 48 353 43 46 01 15  051  002 -001 0.073
%  8323 21 41 456 45 47 45 31 21 48 156 45 47 18 13  056  002 -003 0.073
%  0021  0 03 502 43 58 35 54  0 09 106 43 59 08 59  023  005 -002 0.072
%  0100  0 21 171 45 44 14 05  0 26 121 45 43 40 48  039  001  000 0.072
%  0781 02 34 436 45 12 17 48 02 39 337 45 11 52 20  048  001 -002 0.072
%  0818 02 40 261 45 18 59 45 02 45 061 45 18 34 21  045  003  000 0.072
%  3759 09 24 043 45 02 19 55 09 29 088 45 02 46 08  046  001 -000 0.072
%  4277 10 53 520 43 40 57 52 10 59 279 43 40 25 49  050 -003  001 0.072
%  5404 14 21 475 43 52 18 47 14 25 117 43 51 51 03  040 -002 -004 0.072
%  6998 18 32 556 45 21 08 04 18 38 533 45 21 03 07  059 -001 -002 0.072
%  7377 19 20 274 43 02 54 55 19 25 298 43 03 06 53  034  003  001 0.072
%  1729 05 12 063 43 40  0 37 05 19 084 43 40 05 57  047  005 -007 0.070
%  4623 12 03 152 45 24 10 16 12 08 247 45 24 43 44  040  001 -000 0.070
%  7898 20 34 147 45 24 08 22 20 40 116 45 23 46 26  064  005  005 0.070
%  0370 01 10 403 45 46 03 37 01 15 111 45 45 31 54  050  007  002 0.069
%  0810 02 39 092 45 51 13 54 02 42 334 45 50 48 01  054  003  002 0.069
%  2612 06 53 429 45 35 22 27 06 57 175 45 35 30 28  062 -000  000 0.069
%  3579 08 54 090 43 42 10 44 09  0 383 43 41 46 58  040 -004 -003 0.069
%  5534 14 45 474 43 24 19 29 14 50 158 43 23 54 42  058  002  000 0.069
%  6573 17 33 574 43 61 57 04 17 34 594 43 61 52 30  052  002 -005 0.069
%  2483 06 39 319 43 43 40 37 06 46 443 43 43 34 39  053 -000  002 0.068
%  3775 09 26 103 43 52 08  0 09 32 513 43 51 40 38  032 -010 -005 0.068
%  5447 14 30 195 43 30 10 46 14 34 407 43 29 44 42  045  002  001 0.068
%  8162 21 16 115 43 62 09 43 21 18 347 43 62 35 08  024  001  001 0.068
%  8974 23 35 143 43 77 04 27 23 39 208 43 77 37 57  032 -001  002 0.068
%  0008  0 01 252 43 28 28 11  0 06 367 43 29 01 17  061  004 -002 0.067
%  0033  0 06 105 45 16 01 01  0 11 158 45 15 28 05  049 -001 -003 0.067
%  1355 04 14 454 45 59 32 33 04 16 289 45 59 18 07  044 -001 -002 0.067
%  1780 05 18 352 43 17 17 26 05 24 254 43 17 23  0  050  002 -000 0.067
%  2890 07 28 130 43 32 06 27 07 34 359 43 31 53 18  016 -002 -001 0.067
%  2891 07 28 130 43 32 06 27 07 34 359 43 31 53 18  016 -002 -001 0.067
%  3262 08 13 594 43 27 32 30 08 20 038 43 27 13 03  051 -000 -004 0.067
%  6556 17 30 175 43 12 37 58 17 34 560 43 12 33 36  021  001 -002 0.067
%  7783 20 16 323 43 66 31 56 20 17 311 43 66 51 14  059  005  003 0.067
%  8635 22 36 393 45 47 43 29 22 42 368 45 47 12 39  060   00 -003 0.067
%  0025  0 04 202 45 46 17 57  0 09 246 45 45 44 51  039  001 -002 0.066
%  0244  0 47 058 43 60 34 33  0 53 041 43 61 07 27  048 -001  002 0.066
%  0458 01 30 555 43 40 54 19 01 36 478 43 41 24 20  041 -002 -004 0.066
%  2085 05 51 510 45 14 11 09 05 56 242 45 14 10 04  037 -000  001 0.066
%  5699 15 15 034 45 47 56 57 15 21 481 45 48 19 03  057 -016 -003 0.066
%  0235  0 45 071 45 11 10 58  0 50 075 45 10 38 40  052 -002 -002 0.065
%  1747 05 14 235 45 18 14 12 05 18 504 45 18 07 49  060  004  001 0.065
%  3786 09 26 456 45 40 01 44 09 30 419 45 40 28  0  036 -002  001 0.065
%  5288 14  0 477 45 35 52 41 14 06 408 45 36 22 12  021 -005 -005 0.065
%  8314 21 39 420 43 14 18 58 21 44 315 43 14 46 22  059  003 -001 0.065
%  0652 02 08 302 45 31 11 34 02 12 544 45 30 43 26  053   00  000 0.064
%  1101 03 31 461 43  0 05 04 03 36 523 43  0 24 06  043 -002 -005 0.064
%  1504 04 38 392 45 59 08 25 04 40 183 45 58 56 37  065  001  002 0.064
%  3018 07 41 511 45 33 58 32 07 45 348 45 34 10 23  054 -003  017 0.064
%  3650 09 06 492 43 15 23 57 09 12 175 43 14 59 46  065 -005  002 0.064
%  3684 09 11 449 45 36 59 46 09 15 451 45 37 24 48  046   00 -000 0.064
%  4525 11 42 156 45 29 43 28 11 47 156 45 30 17 13  065 -003 -002 0.064
%  5603 14 58 129 45 24 53 20 15 04 041 45 25 16 55  033 -001 -000 0.064
%  5634 15 02 545 43 25 15 31 15 07 180 43 24 52 09  049  002 -002 0.064
%  1262 03 59 249 43 21 44 21 04 05 202 43 22  0 32  059  002 -001 0.063
%  3430 08 34 452 45 22 19 18 08 39 079 45 22 39 43  050 -002  004 0.063
%  4845 12 40 155 43 39 49 20 12 44 594 43 39 16 44  060 -004  001 0.063
%  7232 19  0 066 45 37 57 10 19 06 523 45 37 48 37  062 -002 -003 0.063
% 
% 080813
%Nearby stars on constellation outlines	ly
% 7377 delta Aquilae			50
% 7557 alpha Aquilae "Altair"		17
% 5235 eta Bootis "Mufrid"		37
% 5340 alpha Bootis "Arcturus"		37
% 2491 alpha Canis Majoris "Sirius"	 9
% 2943 alpha Canis Minoris "Procyon"	11
% 5459 alpha Centauri "Rigel Kentaurus"	 4
% 7957 eta Cephei			47
% 8162 alpha Cephei			49
% 8974 gamma Cephei			45
% 2891 alpha Geminorum "Castor"		52
% 2990 beta Geminorum "Pollux"		34
% 6212 zeta Herculis			35
% 6623 mu Herculis			27
% 7001 alpha Lyrae "Vega"		25
%		55 stars
% 080813
%Filling in the constellations: approximate distances of vertices for planar map
% (Page references \cite{Pasachoff:starPlanets}		55 stars
% Aquila (p.299)	6: alpha, gamma, delta, zeta, theta, lambda
%  alpha - gamma (19 46.2, 10 36, 17)
%  gamma - zeta  (19 05.4, 13 52, 33)
%  alpha - theta (10 11.4, -0 50, 33)
%  zeta - delta
%  theta - delta
%  delta - lambda (19 06.3, -4 53, 96)
% Bootes (p.259)	5: alpha, beta, gamma, delta, eta
%  alpha - delta (15 15.5, 33 19, 37)
%  delta - beta  (15 01.8, 40 24, 37)
%  beta - gamma  (14 32.1, 38 19, 37)
%  gamma - alpha
%  alpha - eta
% Canis Major (p.309)	5: alpha, beta, delta, epsilon, eta
%  alpha - beta  (06 22.7, -17 57, 9)
%  alpha - delta (07 10, -26, 9)	delta is all interpolated by THM
%  delta - epsilon (06 58.6, -28 58, 9)
%  delta - eta   (07 24.1, -29 18, 9)
% Canis Minor (p.285)	2: alpha, beta
%  alpha - beta  (07 27.2, 08 18, 11)
% Centaurus (pp. 315,341,3)
%  9: alpha, beta, gamma, epsilon, zeta, eta, theta, mu, nu
%  alpha - epsilon (13 39.9, -53 28, 4)
%  epsilon - beta  (14 03.8, -60 22, 4)
%  epsilon - gamma (12 41.5, -48 57, 4)
%  epsilon - zeta  (13 55.6, -47 17, 4)
%  zeta - mu       (13 49.6, -42 28, 4)
%  mu - nu         (13 49.5, -41 41, 4)
%  nu - eta        (14 35.5, -42 09, 4)
%  nu - theta      (14 06.7, -36 22, 4)
% Cepheus (p.243)	6: alpha, beta, gamma, zeta, eta iota
%  alpha - eta
%  alpha - beta  (21 28.6, 70 33, 47)
%  beta - gamma
%  gamma - iota  (22 49.6, 66 12, 46)
%  iota - zeta   (22 10.8, 58 12, 47)
%  zeta - alpha
% Gemini (p.251)	 alpha, beta, gamma, delta, epsilon, eta, mu
%  alpha - beta
%  alpha - delta (07 20.1, 21 59, 34)
%  delta - gamma (06 37.7, 16 25, 34)
%  gamma - mu    (06 22.9, 22 31, 52)
%  mu - eta      (06 14.8, 22 31, 52)
%  mu - epsilon  (06 44.0, 25 08, 52)
%  epsilon - alpha
% Hercules (p.261)	11: alpha, beta, delta, epsilon, zeta, eta, theta, iota,
%			    mu, pi, tau
%  mu - delta    (17 15.1, 24 50, 31)
%  delta - epsilon (17, 30 48, 33)	epsilon is all interpolated by THM
%  epsilon - zeta
%  zeta - eta    (16 42.8, 38 56, 35)
%  eta - pi      (17 15.0, 36 48, 32)
%  pi - epsilon
%  zeta - beta   (16 30.3, 21 30, 36)
%  beta - alpha  (17 14,7, 14 23, 29)
%  eta - tau     (16 20, 46, 38)	tau is all interpolated by THM
%  pi - theta    (17 56, 37 30, 25)	theta is all interpolated by THM
%  theta - iota  (17 45, 46, 27)	iota is all interpolated by THM
% Lyra (p.265)		5: alpha, beta, gamma, delta, zeta
%  alpha - zeta  (18 46, 37 40, 25)	zeta is all interpolated by THM
%  zeta - beta   (18 50.0, 33 22, 25)
%  beta - gamma  (18 58.8, 32 41, 25)
%  gamma - delta (18 54, 37, 25)	delta is all interpolated by THM
%  delta - zeta

% drawConst3D % THM 080814	in file drawConst3D.m
constel = cell(56,10);
% RA, DEC, parx for close stars (<= 52 ly) from Yale Bright Star Calatog. Other:
% RA, DEC from \cite{Pasachoff:starPlanets}; parx interpolated from 1--3 close
%		  constel  star      h  m  s   pm  d  m  s  parx
constel( 1,:) = {'Aquila' 'alpha'   19 50 46.9 +1 08 52 06 0.202};
constel( 2,:) = {'Aquila' 'gamma'   19 46 12.0 +1 10 36 00 0.202};
constel( 3,:) = {'Aquila' 'delta'   19 25 29.8 +1 03 06 53 0.072};
constel( 4,:) = {'Aquila' 'zeta'    19 05 24.0 +1 13 52 00 0.137};
constel( 5,:) = {'Aquila' 'theta'   20 11 24.0 -1 00 50 00 0.137};
constel( 6,:) = {'Aquila' 'lambda'  19 06 18.0 -1 04 53 00 0.063};
constel( 7,:) = {'Bootes' 'alpha'   14 15 39.6 +1 19 10 57 0.097};
constel( 8,:) = {'Bootes' 'beta'    15 01 48.0 +1 40 24 00 0.097};
constel( 9,:) = {'Bootes' 'gamma'   14 32 06.0 +1 38 19 00 0.097};
constel(10,:) = {'Bootes' 'delta'   15 15 30.0 +1 33 19 00 0.097};
constel(11,:) = {'Bootes' 'eta'     13 54 41.0 +1 18 23 52 0.108};
constel(12,:) = {'CMajor' 'alpha'   06 45 08.9 -1 16 42 58 0.378};
constel(13,:) = {'CMajor' 'beta'    06 22 42.0 -1 17 57 00 0.378};
constel(14,:) = {'CMajor' 'delta'   07 10 00.0 -1 26 00 00 0.378};
constel(15,:) = {'CMajor' 'epsilon' 06 59 36.0 -1 28 58 00 0.378};
constel(16,:) = {'CMajor' 'eta'     07 27 12.0 -1 29 18 00 0.378};
constel(17,:) = {'CMinor' 'alpha'   07 39 18.1 +1 05 13 30 0.292};
constel(18,:) = {'CMinor' 'beta'    07 27 12.0 +1 08 18 00 0.292};
constel(19,:) = {'Centau' 'alpha'   14 39 36.2 -1 60 50 07 0.750};
constel(20,:) = {'Centau' 'beta'    14 03 48.0 -1 60 22 00 0.750};
constel(21,:) = {'Centau' 'gamma'   12 41 30.0 -1 48 57 00 0.750};
constel(22,:) = {'Centau' 'epsilon' 13 39 54.0 -1 53 28 00 0.750};
constel(23,:) = {'Centau' 'zeta'    13 55 36.0 -1 47 17 00 0.750};
constel(24,:) = {'Centau' 'eta'     14 35 30.0 -1 40 09 00 0.750};
constel(25,:) = {'Centau' 'theta'   14 06 42.0 -1 36 22 00 0.750};
constel(26,:) = {'Centau' 'mu'      13 49 36.0 -1 42 28 00 0.750};
constel(27,:) = {'Centau' 'nu'      13 49 30.0 -1 41 41 00 0.750};
constel(28,:) = {'Cepheu' 'alpha'   21 18 34.7 +1 62 35 08 0.068};
constel(29,:) = {'Cepheu' 'beta'    21 28 36.0 +1 70 33 00 0.068};
constel(30,:) = {'Cepheu' 'gamma'   23 39 20.8 +1 77 37 00 0.068};
constel(31,:) = {'Cepheu' 'zeta'    22 10 48.0 +1 58 12 00 0.068};
constel(32,:) = {'Cepheu' 'eta'     20 45 17.3 +1 61 50 20 0.075};
constel(33,:) = {'Cepheu' 'iota'    22 49 36.0 +1 66 12 00 0.068};
constel(34,:) = {'Gemini' 'alpha'   07 34 35.9 +1 31 53 18 0.067};
constel(35,:) = {'Gemini' 'beta'    07 45 18.9 +1 28 01 34 0.094};
constel(36,:) = {'Gemini' 'gamma'   06 37 42.0 +1 16 25 00 0.094};
constel(37,:) = {'Gemini' 'delta'   07 20 06.0 +1 21 59 00 0.094};
constel(38,:) = {'Gemini' 'epsilon' 06 44 00.0 +1 25 08 00 0.067};
constel(39,:) = {'Gemini' 'eta'     06 14 48.0 +1 22 31 00 0.067};
constel(40,:) = {'Gemini' 'mu'      06 22 54.0 +1 22 31 00 0.067};
constel(41,:) = {'Hercul' 'alpha'   17 14 42.0 +1 14 23 00 0.125};
constel(42,:) = {'Hercul' 'beta'    16 30 18.0 +1 21 30 00 0.099};
constel(43,:) = {'Hercul' 'delta'   17 15 06.0 +1 24 50 00 0.116};
constel(44,:) = {'Hercul' 'epsilon' 17 00 00.0 +1 30 48 00 0.109};
constel(45,:) = {'Hercul' 'zeta'    16 41 17.1 +1 31 36 10 0.102};
constel(46,:) = {'Hercul' 'eta'     16 42 48.0 +1 38 56 00 0.102};
constel(47,:) = {'Hercul' 'theta'   17 56 00.0 +1 37 30 00 0.148};
constel(48,:) = {'Hercul' 'iota'    17 45 00.0 +1 46 00 00 0.136};
constel(49,:) = {'Hercul' 'mu'      17 46 27.5 +1 27 43 15 0.133};
constel(50,:) = {'Hercul' 'pi'      17 15 00.0 +1 36 48 00 0.112};
constel(51,:) = {'Hercul' 'tau'     16 20 00.0 +1 46 00 00 0.093};
constel(52,:) = {'Lyra'   'alpha'   18 36 56.2 +1 38 47 01 0.133};
constel(53,:) = {'Lyra'   'beta'    18 50 00.0 +1 33 22 00 0.133};
constel(54,:) = {'Lyra'   'gamma'   18 58 48.0 +1 32 41 00 0.133};
constel(55,:) = {'Lyra'   'delta'   18 54 00.0 +1 37 00 00 0.133};
constel(56,:) = {'Lyra'   'zeta'    18 46 00.0 +1 37 40 00 0.133};
siz = size(constel);
for k = 1:siz(1)
  RA = pi/180*(360/24*(constel{k,3} + 1/60*(constel{k,4} + 1/60*constel{k,5})));
									% radian
   sD = constel{k,6};							% sign
   %sD = (constel{k,6}=='+')*1 + (constel{k,6}=='-')*(-1);		% sign
   % OR sD = 44 - constel{k,6}
   DEC = pi/180*sD*(constel{k,7} + 1/60*(constel{k,8} + 1/60*constel{k,9}));
									% radian
   ly = 3.26/constel{k,10};						% dist
   xm = ly/105;							% exameters
   X(k) = xm*cos(DEC)*cos(RA);
   Y(k) = xm*cos(DEC)*sin(RA);
   Z(k) = xm*sin(DEC);
end
U = zeros(size(X));
V = zeros(size(Y));
W = zeros(size(Z));
% Aquila alpha - gamma, gamma - zeta, alpha - theta, zeta - delta, theta - delta
%	 delta - lambda		alpha, gamma, delta, zeta, theta, lambda: 123456
U(1) = X(2) - X(1); U(2) = X(4) - X(2); U(4) = X(3) - X(4); U(3) = X(5) - X(3);
U(5) = X(1) - X(5); U(6) = X(3) - X(6);
V(1) = Y(2) - Y(1); V(2) = Y(4) - Y(2); V(4) = Y(3) - Y(4); V(3) = Y(5) - Y(3);
V(5) = Y(1) - Y(5); V(6) = Y(3) - Y(6);
W(1) = Z(2) - Z(1); W(2) = Z(4) - Z(2); W(4) = Z(3) - Z(4); W(3) = Z(5) - Z(3);
W(5) = Z(1) - Z(5); W(6) = Z(3) - Z(6);
% Bootes alpha - delta, delta - beta, beta - gamma, gamma - alpha, alpha - eta
%	alpha, beta, gamma, delta, eta: 7,8,9,10,11
U( 7) = X(10) - X( 7); U(10) = X( 8) - X(10); U( 8) = X( 9) - X( 8);
U( 9) = X( 7) - X( 9); U(11) = X( 7) - X(11);
V( 7) = Y(10) - Y( 7); V(10) = Y( 8) - Y(10); V( 8) = Y( 9) - Y( 8);
V( 9) = Y( 7) - Y( 9); V(11) = Y( 7) - Y(11);
W( 7) = Z(10) - Z( 7); W(10) = Z( 8) - Z(10); W( 8) = Z( 9) - Z( 8);
W( 9) = Z( 7) - Z( 9); W(11) = Z( 7) - Z(11);
% Canis Major alpha - beta, alpha - delta, delta - epsilon, delta - eta
%	alpha, beta, delta, epsilon, eta: 12,13,14,15,16
U(13) = X(12) - X(13); U(14) = X(12) - X(14); U(15) = X(14) - X(15);
U(16) = X(14) - X(16);
V(13) = Y(12) - Y(13); V(14) = Y(12) - Y(14); V(15) = Y(14) - Y(15);
V(16) = Y(14) - Y(16);
W(13) = Z(12) - Z(13); W(14) = Z(12) - Z(14); W(15) = Z(14) - Z(15);
W(16) = Z(14) - Z(16);
% Canis Minor alpha - beta	alpha, beta: 17, 18
U(17) = X(18) - X(17);
V(17) = Y(18) - Y(17);
W(17) = Z(18) - Z(17);
% Centaurus alpha - epsilon, epsilon - beta, epsilon - gamma, epsilon - zeta,
%	    zeta - mu, mu - nu, nu - eta, nu - theta
%	alpha, beta, gamma, epsilon, zeta, eta, theta, mu, nu:
%	 19,    20,   21,    22,      23,   24,  25,   26, 27
U(19) = X(22) - X(19); U(20) = X(22) - X(20); U(21) = X(22) - X(21);
U(22) = X(23) - X(22); U(23) = X(26) - X(23); U(24) = X(27) - X(24);
U(25) = X(27) - X(25); U(26) = X(27) - X(26);
V(19) = Y(22) - Y(19); V(20) = Y(22) - Y(20); V(21) = Y(22) - Y(21);
V(22) = Y(23) - Y(22); V(23) = Y(26) - Y(23); V(24) = Y(27) - Y(24);
V(25) = Y(27) - Y(25); V(26) = Y(27) - Y(26);
W(19) = Z(22) - Z(19); W(20) = Z(22) - Z(20); W(21) = Z(22) - Z(21);
W(22) = Z(23) - Z(22); W(23) = Z(26) - Z(23); W(24) = Z(27) - Z(24);
W(25) = Z(27) - Z(25); W(26) = Z(27) - Z(26);
% Cepheus alpha - eta, alpha - beta, beta - gamma, gamma - iota, iota - zeta,
%	  zeta - alpha
%	alpha, beta, gamma, zeta, eta iota: 28,29,30,31,32,33
U(28) = X(29) - X(28); U(29) = X(30) - X(29); U(30) = X(33) - X(30);
U(33) = X(31) - X(33); U(31) = X(28) - X(31); U(32) = X(28) - X(32);
V(28) = Y(29) - Y(28); V(29) = Y(30) - Y(29); V(30) = Y(33) - Y(30);
V(33) = Y(31) - Y(33); V(31) = Y(28) - Y(31); V(32) = Y(28) - Y(32);
W(28) = Z(29) - Z(28); W(29) = Z(30) - Z(29); W(30) = Z(33) - Z(30);
W(33) = Z(31) - Z(33); W(31) = Z(28) - Z(31); W(32) = Z(28) - Z(32);
% Gemini alpha - beta, beta - delta, delta - gamma, gamma - mu, mu - eta,
%	 mu - epsilon, epsilon - alpha
%	alpha, beta, gamma, delta, epsilon, eta, mu: 34,35,36,37,38,39,40
U(34) = X(35) - X(34); U(35) = X(37) - X(35); U(37) = X(36) - X(37);
U(36) = X(40) - X(36); U(40) = X(38) - X(40); U(38) = X(34) - X(38);
U(39) = X(40) - X(39);
V(34) = Y(35) - Y(34); V(35) = Y(37) - Y(35); V(37) = Y(36) - Y(37);
V(36) = Y(40) - Y(36); V(40) = Y(38) - Y(40); V(38) = Y(34) - Y(38);
V(39) = Y(40) - Y(39);
W(34) = Z(35) - Z(34); W(35) = Z(37) - Z(35); W(37) = Z(36) - Z(37);
W(36) = Z(40) - Z(36); W(40) = Z(38) - Z(40); W(38) = Z(34) - Z(38);
W(39) = Z(40) - Z(39);
% Hercules mu - delta, delta - epsilon, epsilon - zeta, zeta - eta, eta - pi,
%	   pi - epsilon, zeta - beta, beta - alpha, eta - tau, pi - theta,
%	   theta - iota
%	alpha, beta, delta, epsilon, zeta, eta, theta, iota, mu, pi, tau:
%	 41,    42,   43,    44,      45,  46,   47,    48,  49, 50,  51
U(41) = X(42) - X(41); U(42) = X(45) - X(42); U(45) = X(46) - X(45);
U(46) = X(50) - X(46); U(50) = X(44) - X(50); U(44) = X(45) - X(44);
U(49) = X(43) - X(49); U(43) = X(44) - X(43); U(48) = X(47) - X(48);
U(47) = X(50) - X(47); U(51) = X(46) - X(51);
V(41) = Y(42) - Y(41); V(42) = Y(45) - Y(42); V(45) = Y(46) - Y(45);
V(46) = Y(50) - Y(46); V(50) = Y(44) - Y(50); V(44) = Y(45) - Y(44);
V(49) = Y(43) - Y(49); V(43) = Y(44) - Y(43); V(48) = Y(47) - Y(48);
V(47) = Y(50) - Y(47); V(51) = Y(46) - Y(51);
W(41) = Z(42) - Z(41); W(42) = Z(45) - Z(42); W(45) = Z(46) - Z(45);
W(46) = Z(50) - Z(46); W(50) = Z(44) - Z(50); W(44) = Z(45) - Z(44);
W(49) = Z(43) - Z(49); W(43) = Z(44) - Z(43); W(48) = Z(47) - Z(48);
W(47) = Z(50) - Z(47); W(51) = Z(46) - Z(51);
% Lyra alpha - zeta, zeta - beta, beta - gamma, gamma - delta, delta - zeta
%	alpha, beta, gamma, delta, zeta: 52, 53, 54, 55, 56
U(52) = X(56) - X(52); U(56) = X(53) - X(56); U(53) = X(54) - X(53);
U(54) = X(55) - X(54); U(55) = X(56) - X(55);
V(52) = Y(56) - Y(52); V(56) = Y(53) - Y(56); V(53) = Y(54) - Y(53);
V(54) = Y(55) - Y(54); V(55) = Y(56) - Y(55);
W(52) = Z(56) - Z(52); W(56) = Z(53) - Z(56); W(53) = Z(54) - Z(53);
W(54) = Z(55) - Z(54); W(55) = Z(56) - Z(55);
quiver3(X,Y,Z,U,V,W,0,'b.')
xlabel('x'), ylabel('y'), zlabel('z')
hold off