%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialize constants   %
%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all;
UENERGY_IN_JOULE = 1.602192E-19;
ULENGTH_IN_METER = 1.E-10;
UMASS_IN_KG = 1.66053E-27;
UTIME_IN_SECOND = ULENGTH_IN_METER*sqrt(UMASS_IN_KG/UENERGY_IN_JOULE);
UFREQ_IN_HZ = 1./UTIME_IN_SECOND/2/pi;
UFREQ_IN_THZ = UFREQ_IN_HZ / 1E12;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialize Si perfect crystal LDOS   %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
load dos_si.out;
freq = 1:dos_si(2);
Si = sum(dos_si(4:3+dos_si(2)))
clf;
plot((freq-0.5)*dos_si(3)*UFREQ_IN_THZ,dos_si(4:3+dos_si(2))/Si/dos_si(3)/UFREQ_IN_THZ)
;
xlabel ('THz');
ylabel ('1/THz');
hold on;

%%%%%%%%%%%%%%%%
% 4H-SiC (001) %
%%%%%%%%%%%%%%%%
% 896 atoms, 1000 supercell k
load dos4h_Si_001.out;
plot(dos4h_Si_001(:,1),filtfilt([1],[1],dos4h_Si_001(:,2)),'r-.');
axis ([0 35 0 0.18]);

%%%%%%%%%%%%%%%%
% 4H-SiC (100) %
%%%%%%%%%%%%%%%%
% 896 atoms, 1000 supercell k
load dos4h_Si_100.out;
plot(dos4h_Si_100(:,1),filtfilt([1],[1],dos4h_Si_100(:,2)),'g:');
axis ([0 35 0 0.18]);

title ('<001> (dash-dot), <100> (dot) Si LDOS in 4H-SiC compared to 3C-SiC (solid)')

PRINT;
