# 5.8: Numerical Experiment (Star Field)

$$\newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} }$$

$$\newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}}$$

$$\newcommand{\id}{\mathrm{id}}$$ $$\newcommand{\Span}{\mathrm{span}}$$

( \newcommand{\kernel}{\mathrm{null}\,}\) $$\newcommand{\range}{\mathrm{range}\,}$$

$$\newcommand{\RealPart}{\mathrm{Re}}$$ $$\newcommand{\ImaginaryPart}{\mathrm{Im}}$$

$$\newcommand{\Argument}{\mathrm{Arg}}$$ $$\newcommand{\norm}[1]{\| #1 \|}$$

$$\newcommand{\inner}[2]{\langle #1, #2 \rangle}$$

$$\newcommand{\Span}{\mathrm{span}}$$

$$\newcommand{\id}{\mathrm{id}}$$

$$\newcommand{\Span}{\mathrm{span}}$$

$$\newcommand{\kernel}{\mathrm{null}\,}$$

$$\newcommand{\range}{\mathrm{range}\,}$$

$$\newcommand{\RealPart}{\mathrm{Re}}$$

$$\newcommand{\ImaginaryPart}{\mathrm{Im}}$$

$$\newcommand{\Argument}{\mathrm{Arg}}$$

$$\newcommand{\norm}[1]{\| #1 \|}$$

$$\newcommand{\inner}[2]{\langle #1, #2 \rangle}$$

$$\newcommand{\Span}{\mathrm{span}}$$ $$\newcommand{\AA}{\unicode[.8,0]{x212B}}$$

$$\newcommand{\vectorA}[1]{\vec{#1}} % arrow$$

$$\newcommand{\vectorAt}[1]{\vec{\text{#1}}} % arrow$$

$$\newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} }$$

$$\newcommand{\vectorC}[1]{\textbf{#1}}$$

$$\newcommand{\vectorD}[1]{\overrightarrow{#1}}$$

$$\newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}}$$

$$\newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}}$$

$$\newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} }$$

$$\newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}}$$

$$\newcommand{\avec}{\mathbf a}$$ $$\newcommand{\bvec}{\mathbf b}$$ $$\newcommand{\cvec}{\mathbf c}$$ $$\newcommand{\dvec}{\mathbf d}$$ $$\newcommand{\dtil}{\widetilde{\mathbf d}}$$ $$\newcommand{\evec}{\mathbf e}$$ $$\newcommand{\fvec}{\mathbf f}$$ $$\newcommand{\nvec}{\mathbf n}$$ $$\newcommand{\pvec}{\mathbf p}$$ $$\newcommand{\qvec}{\mathbf q}$$ $$\newcommand{\svec}{\mathbf s}$$ $$\newcommand{\tvec}{\mathbf t}$$ $$\newcommand{\uvec}{\mathbf u}$$ $$\newcommand{\vvec}{\mathbf v}$$ $$\newcommand{\wvec}{\mathbf w}$$ $$\newcommand{\xvec}{\mathbf x}$$ $$\newcommand{\yvec}{\mathbf y}$$ $$\newcommand{\zvec}{\mathbf z}$$ $$\newcommand{\rvec}{\mathbf r}$$ $$\newcommand{\mvec}{\mathbf m}$$ $$\newcommand{\zerovec}{\mathbf 0}$$ $$\newcommand{\onevec}{\mathbf 1}$$ $$\newcommand{\real}{\mathbb R}$$ $$\newcommand{\twovec}[2]{\left[\begin{array}{r}#1 \\ #2 \end{array}\right]}$$ $$\newcommand{\ctwovec}[2]{\left[\begin{array}{c}#1 \\ #2 \end{array}\right]}$$ $$\newcommand{\threevec}[3]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \end{array}\right]}$$ $$\newcommand{\cthreevec}[3]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \end{array}\right]}$$ $$\newcommand{\fourvec}[4]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}$$ $$\newcommand{\cfourvec}[4]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}$$ $$\newcommand{\fivevec}[5]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}$$ $$\newcommand{\cfivevec}[5]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}$$ $$\newcommand{\mattwo}[4]{\left[\begin{array}{rr}#1 \amp #2 \\ #3 \amp #4 \\ \end{array}\right]}$$ $$\newcommand{\laspan}[1]{\text{Span}\{#1\}}$$ $$\newcommand{\bcal}{\cal B}$$ $$\newcommand{\ccal}{\cal C}$$ $$\newcommand{\scal}{\cal S}$$ $$\newcommand{\wcal}{\cal W}$$ $$\newcommand{\ecal}{\cal E}$$ $$\newcommand{\coords}[2]{\left\{#1\right\}_{#2}}$$ $$\newcommand{\gray}[1]{\color{gray}{#1}}$$ $$\newcommand{\lgray}[1]{\color{lightgray}{#1}}$$ $$\newcommand{\rank}{\operatorname{rank}}$$ $$\newcommand{\row}{\text{Row}}$$ $$\newcommand{\col}{\text{Col}}$$ $$\renewcommand{\row}{\text{Row}}$$ $$\newcommand{\nul}{\text{Nul}}$$ $$\newcommand{\var}{\text{Var}}$$ $$\newcommand{\corr}{\text{corr}}$$ $$\newcommand{\len}[1]{\left|#1\right|}$$ $$\newcommand{\bbar}{\overline{\bvec}}$$ $$\newcommand{\bhat}{\widehat{\bvec}}$$ $$\newcommand{\bperp}{\bvec^\perp}$$ $$\newcommand{\xhat}{\widehat{\xvec}}$$ $$\newcommand{\vhat}{\widehat{\vvec}}$$ $$\newcommand{\uhat}{\widehat{\uvec}}$$ $$\newcommand{\what}{\widehat{\wvec}}$$ $$\newcommand{\Sighat}{\widehat{\Sigma}}$$ $$\newcommand{\lt}{<}$$ $$\newcommand{\gt}{>}$$ $$\newcommand{\amp}{&}$$ $$\definecolor{fillinmathshade}{gray}{0.9}$$

With Earth at the origin, we may specify the star positions for the Big and Little Dippers in three-dimensional homogeneous coordinates. With light years as our units, we have G$$^T$$ =

1.5441      -1.2064     153.0875    1.0000
-1.0386     8.4588      142.7458    1.0000
-8.7955     26.2870     198.0698    1.0000
-12.8807    19.0964     106.5383    1.0000
-18.8926    17.461      90.6185     1.0000
-45.1364    54.2706     215.1148    1.0000
-9.6222     20.1734     88.0062     1.0000
-33.7097    -8.4048     64.6574     1.0000
-33.7144    -8.7191     52.3806     1.0000
-43.8531    -1.1483     59.7193     1.0000
-36.1086    2.3667      55.7927     1.0000
-34.2309    8.2181      52.1260     1.0000
-30.7876    11.8182     46.9810     1.0000
-61.8581    31.5183     80.7615     1.0000

To make use of this data, we need a function to display it. Enter and save the following generalization of the function from Demo 1 in "Vector Graphics: Two Dimensional Image Representation". Call it vhgraph.m.

function vhgraph (P,L,A,ps,ls);
% vhgraph(P,L,A,PS,LS) graphs images whose points are
% stored in P and whose lines are stored in L. The points
% in P must be in homogeneous coordinates in either 2 or 3
% dimensions, with each column of P representing a point.
% The lines are coded in L, with each column of L containing
% 2 integers pointing to a pair of points in P to be
% connected by a line segment. If A is present, the points
% in P are transformed to A*P beforre graphing. For 3D data
% points, only the first two coordinates are graphed, so A
% should include the desired projection from 3D to 2D. The
% point symbol may be specified in PS and the line type in
% LS, if desired
%
% Richard T. Behrens, October 1989.
%
% The first dection of the program determines the sizes of
% all the input matrices and chacks that they make sense.
[mp,np] = size(P);
if (nargin > 1)
[ml,nl] = size(L);
else
ml = 2; nl = 0;
end
if (nargin <=2)
A = eye(mp);
end
[mA,nA] = size(A);
if ((mp~= nA) | (ml~=2))
error('Incompatible sizes of input matrices.")
end
if (nargin <= 3);
ps = '*';
end

P = A*P;    % Performs the transformation A on the points
% (effect is only local to this function).

% The next section contains a loop that renormalizes the
% homogeneous coordinates of the points if necessary.
renorm = find((P(mA,:)~= 1));
if~isempty(renorm)
for i = 1:length(renorm)
P(:,renorm(i)) = P(:,renorm(i))/P(mA,renorm(i));
end
end

% The next program line sets a fixed scale output window
% from -50 to 50 in both x and y directions on the screen.
% For automatic scaling to include all points of the
% image, we could use instead the line q = min(min(P));
% r = max(max(P));
q = -50; r = 50;
axis([q r q r])
axis('square')
plot(P(1,:),P(2,:),ps)  % Plots the points with symbol ps.
hold on                 % Saves the points while we plot
for i = 1:nl            % lines with line type LS.
plot([P(1,L(1,i))P(1,L(2,i))],..
[P(2,L(1,i))P(2,L(2,i))],ls)
end
hold off

Enter the point matrix given at rhe beginning of this section (and take its transpose to put it in the usual form). Also enter the line matrix from Demo 1 in "Vector Graphics: Two Dimensional Image Representation". Save these two matrices and try looking at the image

≫ save dip3d
≫ vhgraph(G,H)

No dippers in sight, right? Without specifying a transformation matrix $$A$$, we have defaulted to looking down on the $$x-y$$ plane from $$z=\infty$$ (a parallel projection). This is how the constellations would look from a distant galaxy (say, a billion light years north of here) through an enormous telescope. We need a perspective view from the origin (Earth), but first we need a set of functions to give us the fundamental operators with which we can build the desired projection.

Take R$$_y(\theta)$$ as an example. The function to build it looks like

 function Ry = vhry (theta) ;
% Rotation matrix about y-axis for 3-D homogeneous
% coordinates.
Ry = eye(4);
Ry(1,1) = cos(theta);
Ry(3,3) = cos(theta);
Ry(3,1) = -sin(theta);
Ry(1,3) = sin(theta);

Enter and save hry.m as given. Then write functions for

$\begin{gathered} \mathrm{R}_{x}(\theta) \quad vhrx . \mathrm{m} \\ \mathrm{R}_{z}(\theta) \quad vhrz. \mathrm{m} \\ \mathrm{S}\left(s_{x}, s_{y}, s_{z}\right) \quad vhs. \mathrm{m} \\ \mathrm{T}\left(t_{x}, t_{y}, t_{z}\right) \quad vht. \mathrm{m} \\ \mathrm{M}(d) \quad vhm. \mathrm{m} \end{gathered} \nonumber$

Useful MATLAB functions for this task include zeros, eye, and diag.

Now build and use a perspective projection with viewpoint Earth and projection plane at $$z=1000$$ behind the dippers:

1. translate Earth to $$z=-1000$$ so that the projection plane coincides with the $$x-y$$ plane: T(O,0,−1000);
2. use the fundamental perspective projection: M(-1000); and
3. translate back: T(O,0,1000)
≫ A = vht(0,0,1000) * vhm(-1000) * vht(0,0,-1000)
≫ vhgraph(G,H,A)

Oops! Now the image is too big; it's mostly off the screen. Scale it down and have another look:

≫ A = vhs(.06,.06,.06) * A
≫ vhgraph(G,H,A)

Now the view should look familiar. Leave $$A$$ as it is now:

≫ save dip3d

Experiment with scale and rotation about the z-axis. For example, try

≫ vhgraph(G,H,vhrz(pi/2) * A)

The two-dimensional star positions given in Demo 1 in "Vector Graphics: Two Dimensional Image Representation" were obtained from the three-dimensional positions with the composite operator A you are now using. To compare the two, type

≫ Gnew = A * G

≫ for i = 1:14
Gnew(;,i) = Gnew(:,i)/Gnew(4,i);
end

≫ Gnew

and compare the $$x$$ and $$y$$ coordinates with those of Demo 1 in "Vector Graphics: Two Dimensional Image Representation".

Astronomers give star positions in equatorial coordinates using right ascension, declination, and distance. The following function converts equatorial coordinates, which are spherical, to Cartesian coordinates with the z-axis pointing north, the x-axis pointing at the vernal (Spring) equinox in the constellation Pisces, and the y-axis pointing toward the Winter solstice in the constellation Opheuchus.

function v = starxyz(rah,ram,decd,decm,dist)
% starxyz is the cartesian coordinates of a star whose
% spherical coordinates ( e.g. from a star catalog) are
%
%          rah     right ascension hours
%          ram     right ascension minutes
%          decd    declination degrees
%          decm    declination minutes (should be negative
%                  if decd is negative)
%          dist    distance (light years)
%
phi = (pi/180) * (decd + decm/60);
theta = (pi/12) * (rah + ram/60);
r = dist;
v = [r * cos(phi) * cos(theta); -r * cos(phi) * sin(theta);
r * sin(phi)];

Have you ever wondered what the constellations would look like from other places in the galaxy? We will soon see the answer. First we will view the dippers from Alpha Centauri, the nearest star, whose coordinates are

 -1.5680     1.3157      -3.6675.

We will look toward the centroid of the fourteen stars in the dippers, located at

-26.3632      12.8709      100.4714.

To get the desired view, we must

1. translate the viewpoint to the origin;
2. rotate the centroid (direction to look) to the z-axis–note that the centroid will have new coordinates after step (1); and
3. apply the composite A=S(.06,.06,.06)*T(0,0,1000)*M(-1000)*T(0,0,-1000) (as used to view from Earth).

Write a function vhz.m based on Exercise 5 from "Vector Graphics: Three-Dimensional Homogeneous Coordinates" to accomplish step (2). Test it on several random points to make sure it works right. Now write a general perspective projection function called vhview.m. The function vhview should accept as inputs two vectors, the first specifying the viewpoint and the second the point to look toward. Its output should be a composite operator that performs all three of the preceding steps.

Now we want to look toward the centroid of the dippers from Alpha Centauri. To do so, enter the vectors for the two points of interest and construct the view like this:

≫ acentauri = [-1.5680; 1.3157; -3.6675]
≫ centroid = [-26.3632; 12.8709; 100.4714]
≫ A = vhview(acentauri,centroid)
≫ vhgraph(G,H,A)
≫ title('From Alpha Centauri')

The farther we move from Earth, the more distorted the dippers will look in general. It should be easy now to view them from any desired location. Just choose a viewpoint, recalculate the composite operator for that viewpoint using vhview, and use vhgraph. Follow this procedure to view the dippers from each of the stars in the following list. You will need to use starxyz first to convert their coordinates.

 Star Right Ascension Declination Distance (ly) Alpha Centauri 14h 40m -60° 50' 4.2 Sirius 6h 45m -16°43' 9.5 Arcturus 14h 16m 19°11' 16.6 Pollux 7h 45m 28°02' 35.9 Betelgeuse 5h 55m 7°24' 313.5

Of course, star viewing is not the only application of vector graphics. Do some experiments with the unit cube (see Exercise 2 from "Vector Graphics: Three-Dimensional Homogeneous Coordinates"). View the cube from location (4,3,2) looking toward the origin using the procedure just outlined for stars. You may need to adjust the scaling to get a meaningful view.

This page titled 5.8: Numerical Experiment (Star Field) is shared under a CC BY 3.0 license and was authored, remixed, and/or curated by Louis Scharf (OpenStax CNX) via source content that was edited to the style and standards of the LibreTexts platform.