function [strobe, reveal, daylight, halogen] =  illuminance
%Psych 221
%Infrared Photography-Final Project
%Bill Peter
%Charles Chang-I Wang
close all
clear all
%This file plots spectral illuminance of different light sources

%Start with Incadescent (Strobe Light)
%Scale by 1260 Lumens, 25% duty cycle
output = .25
wavelength = [400e-9:1e-9:900e-9];
c = 3e8;
h = 6.626e-34
k = 1.38e-23
Temperature = 2700 %in Kelvin

front = ((8*pi*h*c)./(wavelength.^5));

bottom = exp((h*c)./(Temperature*k*wavelength)) -1;

strobe = output.*(front .* (1./bottom));

plot(wavelength, strobe,'linewidth',3)
grid on
title('Spectra for Strobe Light (2700K)','fontsize',14)
xlabel('wavelength','fontsize',14)
ylabel('Spectral Energy (J/m^3)','fontsize',14)


%Now make filter for Reveal glass bulb
%Takes out yellow light so make notch around 550nm
%amount decreased enter into scale
%will have to take 1 - this to get transfer function
%Operates at 1260 lumens
scale(1:150) = 0;
scale(151:160)= linspace(0,.01,10);
scale(161:170)= linspace(.01,.05,10);
scale(171:200)= linspace(.06,.44,30);
scale(201:210)= linspace(.44,.49,10);
scale(211:220)= linspace(.49,.5,10);
scale(221:501) = 0.5;
adjust_scale = 1-scale;

figure
plot(wavelength,adjust_scale)
figure
reveal = (strobe.*adjust_scale)/output;
scale_factor = sum(reveal)/(1e5)
plot(wavelength, reveal,'linewidth',3)
grid on
title('Spectra for Reveal 100 Bulb (2700K)','fontsize',14)
xlabel('wavelength','fontsize',14)
ylabel('Spectral Energy (J/m^3)','fontsize',14)





figure
%Now for Halogen bulb
%Operated at 870 lumens
%must bring scale relative to incadescent output
wavelength = [400e-9:1e-9:900e-9];
c = 3e8;
h = 6.626e-34;
k = 1.38e-23;
Temperature = 4000; %in Kelvin

front = ((8*pi*h*c)./(wavelength.^5));

bottom = exp((h*c)./(Temperature*k*wavelength)) -1;

halogen = (front .* (1./bottom));

scale_halogen = sum(halogen)*(1e-5)
rescale = ((scale_factor)/(scale_halogen))
halogen = rescale*halogen;
plot(wavelength, halogen,'linewidth',3)
grid on
title('Spectra for Halogen Bulb (4000K)','fontsize',14)
xlabel('wavelength','fontsize',14)
ylabel('Spectral Energy (J/m^3)','fontsize',14)

figure
%Now for Daylight
%Modeled as at 2,000 lumens
%Numbers not to be compared directly with others since lumens not known
output = 2;
wavelength = [400e-9:1e-9:900e-9];
c = 3e8;
h = 6.626e-34;
k = 1.38e-23;
Temperature = 6000; %in Kelvin

front = ((8*pi*h*c)./(wavelength.^5));

bottom = exp((h*c)./(Temperature*k*wavelength)) -1;

daylight = output.*(front .* (1./bottom));

plot(wavelength, daylight,'linewidth',3)
grid on
title('Spectra for Daylight (6000K)','fontsize',14)
xlabel('wavelength','fontsize',14)
ylabel('Spectral Energy (J/m^3)','fontsize',14)
