new; 
clf;
 
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;
 
load dos_si.out;
freq = 1:dos_si(2);
Si = sum(dos_si(4:3+dos_si(2)))
 
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');
hold on;

load dos.out;
plot(dos(:,1),dos(:,2),'--');

axis([ 0 40 -0.01 0.15]);
PRINT;

clf;

load dos_c.out;
freq = 1:dos_c(2);
C = sum(dos_c(4:3+dos_c(2)))
 
plot((freq-0.5)*dos_c(3)*UFREQ_IN_THZ,dos_c(4:3+dos_c(2))/C/dos_c(3)/UFREQ_IN_THZ);
xlabel ('THz');
hold on;

load doc.out;
plot(doc(:,1),doc(:,2),'--');
axis ([0 40 0 0.6]);






