#Doping profile

yslice = -50e-9;

npts = 1001;
xmin = getnamed("CHARGE simulation region","x min");
xmax = getnamed("CHARGE simulation region","x max");

x = linspace(xmin,xmax,npts);
y = matrix(npts,1);
y(1:npts) = yslice;

# Get the triangular mesh data

vtx3D = getdata("CHARGE","doping.vertices");
nvtx = length(vtx3D)/3;
vtx = [vtx3D(1:nvtx,1), vtx3D(1:nvtx,3)];
tri = getdata("CHARGE","doping.elements");

# Interpolate the doping profile to rectangular grid

NA = pinch(getdata("CHARGE","doping","NA"));
ND = pinch(getdata("CHARGE","doping","ND"));

NA = pinch(NA(1:nvtx,1));
ND = pinch(ND(1:nvtx,1));

NAi = interptri(tri,vtx,NA,x,yslice);
NDi = interptri(tri,vtx,ND,x,yslice);

plot(x*1e6,log10(abs(NDi - NAi)),"X (um)", "Net Doping (1/cm^3, log scale)");