#######  Script to generate radial doping profile  ####

# Inputs

NA = -1e19;  # donor concentration; unit /cm^3.  Use negative value for acceptors
ND = 1e19;
P = 101;  # number of points along x and y direction

# rectangular grid

x = linspace(-0.2,0.2,P)*1e-6;
y = linspace(-0.2,0.2,P)*1e-6;
z = linspace(-0.53,3.57,2)*1e-6;  # one pont at the bottom and one at the top

N = matrix(P,P,2);

R_outer = 0.095e-6;  # outer radius
R_inner = 0.0251e-6;  # inner radius

for (i=1:P) {
    for (j=1:P) {
        if (sqrt(x(i)^2+y(j)^2)>R_outer) {    #outside
             N(i,j,1:2) = 0;} 
        if (sqrt(x(i)^2+y(j)^2)<R_inner) {    #core
              N(i,j,1:2) = ND*1e6;}
        if ((sqrt(x(i)^2+y(j)^2)<R_outer)and(sqrt(x(i)^2+y(j)^2)>R_inner)) {#shell
              N(i,j,1:2) = NA*1e6;}      
    }
}
# plot profile on XY plane

image(x,y,N(1:P,1:P,1),'x (m)','y (m)','N (/m^3)');

# save data in .mat file

matlabsave("charge",N,x,y,z);