####################################################
# File: CMOS_psf_analysis.lsf
#
# Description:  This file will plot some results from 
#  the frequency sweep object "sweep x" in 
#  CMOS_psf.fsp
#
# Copyright 2011 Lumerical Solutions
####################################################


sname = "sweep x";

# image the incoherent sum of source beams
xs = getsweepdata(sname,"x_source");
ys = getsweepdata(sname,"y_source");
Pz_source = getsweepdata(sname,"Pz_source");
image(xs*1e6,ys*1e6,-0.5*real(Pz_source),"x (microns)","y (microns)","Source power");

# image the Poynting vector at the surface
x = getsweepdata(sname,"x");
y = getsweepdata(sname,"y");
Pz = getsweepdata(sname,"Pz");
image(x*1e6,y*1e6,-0.5*real(Pz),"x (microns)","y (microns)","Power at Si surface");
image(x*1e6,y*1e6,-0.5*real(Pz),"x (microns)","y (microns)","Power at Si surface","logplot");

# get results for each pixel
Power_blue = getsweepdata(sname,"Pblue");
Power_red = getsweepdata(sname,"Pred");
Power_green = getsweepdata(sname,"Pgreen");
Power_all = getsweepdata(sname,"Pall");

# display to screen
?"power absorbed in Si layer: " + num2str(round(Power_all*1000)/10) + "%";
?"power absorbed in all pixels: " + num2str(round((sum(Power_red+Power_green+Power_blue))*1000)/10) + "%";
?"power to red pixels (%): ";
?num2str(round(Power_red*1000)/10);
?"power to green pixels (%): "; 
?num2str(round(Power_green*1000)/10);
?"power to blue pixels (%): "; 
?num2str(round(Power_blue*1000)/10);

# calculate the MTF for green
n = size(Power_green);
x0 = (1:n(1))*4e-6;
y0 = (1:n(2))*4e-6;
sf_x = fftk(x0,3,2^10)/(2*pi);
sf_y = fftk(y0,3,2^10)/(2*pi);
MTF = abs(fft(Power_green,3,2^10));
MTF=MTF/max(MTF);
image(sf_x*1e-3,sf_y*1e-3,MTF,"x cycles per mm","y cycles per mm","MTF");