% Example program to call spharm.m
%
% This program makes two plots. The first is a rectangular plot of the
% selected harmonic (red=positive, blue=negative). The second plot shows
% the function mapped onto a sphere.
%

% Use the n=3, m=2 harmonic

n = 3;
m = 2;

% With Cnm = 1 ; Snm = 0

Cnm = 1;
Snm = 0;

% Define a regular grid of theta (colatitude) and phi (longitude)

theta = linspace(0, pi, 101);
phi = linspace(0, 2*pi, 201);

% Calculate the harmonic function a grid of points defined by the theta and
% phi values

vals = spharm(n, m, Cnm, Snm, theta, phi);

% Now make figures

% first figure: a rectangular plot
%  Note the transformation of the axes.
%  x-axis: longitude given in degrees
%  y-axis: latitude in degrees
%          recall that latitude = 90 - colatitude

figure
pcolor(phi*180/pi, 90-theta*180/pi, vals);
shading('flat');

% second figure: map the harmonic onto a sphere.
%   We need to reverse the order of the rows, so that the northern
%   hemisphere (theta = 0 to pi/2) ends up on the top of the sphere.

[nx ny] = size(vals);
for i=1:nx,
    vals2(i,:) = vals(nx+1-i,:);
end

figure
sphere; h = findobj(gcf, 'Type', 'surface');
set(h, 'CData', vals2, 'FaceColor', 'texturemap')
axis equal;
ylabel('long = 0')
