% This file plots the photon flux densities versus wavelength and frequency
%  for the temperatures of the 3 stars mentioned in the notes.
%  It also plots bands of blue, green, yellow, and red colors. 

% Define vector of wavelengths:
lam = [0.2:0.01:2];   % wavelength [micron]

% Common constant:
c =          % Find online the value of the speed of light and record it to 2 significant figures
             %  in the units of micron per picosecond.
             %  The reason why we will use these units is because we will express 
             %  wavelengths in microns and frequencies in terahertz. 
h =          % Find the value of the Plack constant h (not h-bar) and express it in Joule*picosecond.
k =          % Find the value of the Boltzmann constant and express it in Joule/kiloKelvin.
             
hconk = h*c/k; % h*c/k  [micron*kiloKelvin],  
               % Note: Do not multiply these constants yourself! Matlab will do it better than you. 
             
%  --------------------- Define photon densities as a function of wavelength ------------------------
%  ------------  using Eq. (13.41) in Lecture 13, while omitting the factor const'*T^4.   -----------
%
% Betelgeuse:
TB = 3.4;               %temperature [kiloKelvin]
yB = hconk./(TB*lam);   % nondimensional parameter  y  of eq. (13.41)
NlamB = (yB.^4)./(exp(yB)-1);  % photon spectral density

% Sun:
TSun =               %temperature [kiloKelvin]
% Complete the rest by analogy with Betelgeuse:
ySun = 
NlamSun = 

% Sirius:
TS = 
% complete the rest by analogy with Betelgeuse and Sun


% Photon densities in terms of frequency:
nu = c./lam;          % frequency in THz
honk =                % h/k in ps*kiloKelvin 
%
% Betelgeuse:
zB = (honk/TB)*nu;
NnuB =                % the expression for N(nu), found at the end of Sec. 13.5.
                      % Notes:  1) omega = 2*pi*nu, but hbar = h/(2*pi), so that hbar*omega = h*nu.
                      %         2) Enter only the part of NnuB that depends on zB; omit all other
                      %            multiplicative constants, similarly to how this is done for NlamB. 

% Sun:
zSun = 
% Complete the rest by analogy with Betelgeuse.

% Sirius:
zS = 
% Complete the rest by analogy with Betelgeuse.


%  This is a cosmetic part that will plot bands of colors.
%  Find the colors' wavelengths in Table a link to which is posted next to Lecture 13.
% 
lamBmin =        % minimum wavelength of the Blue band, microns
lamBmax =        % maximum wavelength of the Blue band, microns

lamGmin =        % green
lamGmax = 

lamYmin =        % yellow
lamYmax = 

lamRmin =        % red
lamRmax = 


% Plot photon spectral densities versus wavelength along with color bands:
%  (you do not need to change anything in this part)
figure(131);
plot(lam,NlamB,'r', lam,NlamSun,'y', lam,NlamS,'b','linewidth',3)
axis([min(lam) max(lam) 0 1.2*max(NlamB)])
xlabel('wavelength,  \mum','fontsize',18)
ylabel('N(\lambda,T)','fontsize',18)
patch([lamBmin lamBmin lamBmax lamBmax], [0 1.2*max(NlamB) 1.2*max(NlamB) 0], [84, 176, 239]/255) 
patch([lamGmin lamGmin lamGmax lamGmax], [0 1.2*max(NlamB) 1.2*max(NlamB) 0], [102, 255, 102]/255) 
patch([lamYmin lamYmin lamYmax lamYmax], [0 1.2*max(NlamB) 1.2*max(NlamB) 0], [255, 255, 102]/255) 
patch([lamRmin lamRmin lamRmax lamRmax], [0 1.2*max(NlamB) 1.2*max(NlamB) 0], [255, 90, 100]/255) 

% Plot photon spectral densities versus frequency along with color bands:
%  (you do not need to change anything in this part)
figure(132);
plot(nu,NnuB,'r', nu,NnuSun,'y', nu,NnuS,'b','linewidth',3)
axis([min(nu) 1000 0 1.2*max(NnuB)])
xlabel('frequency,  THz','fontsize',18)
ylabel('N(\nu,T)','fontsize',18)
patch(c./[lamBmin lamBmin lamBmax lamBmax], [0 1.2*max(NnuB) 1.2*max(NnuB) 0], [84, 176, 239]/255) 
patch(c./[lamGmin lamGmin lamGmax lamGmax], [0 1.2*max(NnuB) 1.2*max(NnuB) 0], [102, 255, 102]/255) 
patch(c./[lamYmin lamYmin lamYmax lamYmax], [0 1.2*max(NnuB) 1.2*max(NnuB) 0], [255, 255, 102]/255) 
patch(c./[lamRmin lamRmin lamRmax lamRmax], [0 1.2*max(NnuB) 1.2*max(NnuB) 0], [255, 90, 100]/255) 


