Nanoexplorer

Declarations

Miscellaneous

npix=128; nh=n/2; nq=n/4;

aOff[array_]:=RotateLeft[array,{nh,nh}]

AbsSqOf[z_]:=((Abs[z])^2)

ArrayPwr[Array_]:=Map[AbsSqOf,Array,{2}]

Display Routines

DotsZ[CoordListZ_,xx_,yy_,zz_,srad_]:=Module[{pts},
		dots2[{x_,y_,z_,Z_}]:=Point[{x+xx,y+yy,z+zz}];
	  pts = Table[CoordListZ[[i]],{i,1,Length[CoordListZ]}];
		pts=Map[dots2,pts];
	pts=Table[
{
PointSize[srad*((ElementData[CoordListZ[[i,4]],"AtomicRadius"])/1000)],ColorData["Atoms",CoordListZ[[i,4]]],
Opacity[.5],
pts[[i]]
},{i,1,Length[pts]}];

	{pts}];

Fire5={{0,Black},{.1,Darker[Blue,.35]},{.25,Darker[Magenta,.25]},{.4,Darker[Red,.15]},{.6,Darker[Orange,0.05]},{.8,Yellow},{1,White}};

Specimens

dccxyz4={{-5.43`,0,0,14},{-5.43`,0,5.43`,14},{-5.43`,5.43`,0,14},{-5.43`,5.43`,5.43`,14},{0,-5.43`,0,14},{0,-5.43`,5.43`,14},{0,0,-5.43`,14},{0,0,0,14},{0,0,5.43`,14},{0,0,10.86`,14},{0,5.43`,-5.43`,14},{0,5.43`,0,14},{0,5.43`,5.43`,14},{0,5.43`,10.86`,14},{0,10.86`,0,14},{0,10.86`,5.43`,14},{5.43`,-5.43`,0,14},{5.43`,-5.43`,5.43`,14},{5.43`,0,-5.43`,14},{5.43`,0,0,14},{5.43`,0,5.43`,14},{5.43`,0,10.86`,14},{5.43`,5.43`,-5.43`,14},{5.43`,5.43`,0,14},{5.43`,5.43`,5.43`,14},{5.43`,10.86`,0,14},{10.86`,0,0,14},{10.86`,0,5.43`,14},{10.86`,5.43`,0,14},{-2.715`,-2.715`,0,14},{-2.715`,-2.715`,5.43`,14},{-2.715`,2.715`,-5.43`,14},{-2.715`,2.715`,0,14},{-2.715`,2.715`,5.43`,14},{-2.715`,8.145`,0,14},{-2.715`,8.145`,5.43`,14},{2.715`,-2.715`,-5.43`,14},{2.715`,-2.715`,0,14},{2.715`,-2.715`,5.43`,14},{2.715`,2.715`,-5.43`,14},{2.715`,2.715`,0,14},{2.715`,2.715`,5.43`,14},{2.715`,2.715`,10.86`,14},{2.715`,8.145`,-5.429999`,14},{2.715`,8.145`,0,14},{2.715`,8.145`,5.43`,14},{8.145`,-2.715`,0,14},{8.145`,-2.715`,5.43`,14},{8.145`,2.715`,-5.429999`,14},{8.145`,2.715`,0,14},{8.145`,2.715`,5.43`,14},{8.145`,8.145`,0,14},{8.145`,8.145`,5.4300000000000015`,14},{-5.43`,-2.715`,2.7149994`,14},{-5.43`,2.715`,-2.715`,14},{-5.43`,2.715`,2.7149994`,14},{-5.43`,2.715`,8.145`,14},{-5.429999`,8.145`,2.715`,14},{0,-2.715`,-2.715`,14},{0,-2.715`,2.715`,14},{0,-2.715`,8.145`,14},{0,2.715`,-2.715`,14},{0,2.715`,2.715`,14},{0,2.715`,8.145`,14},{0,8.145`,-2.7149994`,14},{0,8.145`,2.715`,14},{0,8.145`,8.145`,14},{5.43`,-2.715`,-2.7149994`,14},{5.43`,-2.715`,2.715`,14},{5.43`,-2.715`,8.145`,14},{5.43`,2.715`,-2.7149994`,14},{5.43`,2.715`,2.715`,14},{5.43`,2.715`,8.145`,14},{5.43`,8.145`,-2.714999`,14},{5.43`,8.145`,2.715`,14},{5.431`,8.145`,8.145`,14},{10.86`,2.715`,2.715`,14},{-2.715`,-5.43`,2.7149994`,14},{-2.715`,0,-2.715`,14},{-2.715`,0,2.715`,14},{-2.715`,0,8.145`,14},{-2.7149994`,5.43`,-2.7149994`,14},{-2.7149994`,5.43`,2.715`,14},{-2.7149994`,5.43`,8.145`,14},{2.7149994`,-5.43`,-2.715`,14},{2.7149994`,-5.43`,2.7149994`,14},{2.7149994`,-5.43`,8.145`,14},{2.715`,0,-2.715`,14},{2.715`,0,2.715`,14},{2.715`,0,8.145`,14},{2.715`,5.43`,-2.7149994`,14},{2.715`,5.43`,2.715`,14},{2.715`,5.43`,8.145`,14},{2.715`,10.86`,2.715`,14},{8.145`,-5.43`,2.715`,14},{8.145`,0,-2.7149994`,14},{8.145`,0,2.715`,14},{8.145`,0,8.145`,14},{8.145`,5.43`,-2.714999`,14},{8.145`,5.43`,2.715`,14},{8.145`,5.43`,8.145`,14},{-4.0725`,-4.0725`,1.3574995`,14},{-4.0725`,-4.0725`,6.787499`,14},{-4.0725`,1.3575`,-4.0725`,14},{-4.0725`,1.3575`,1.3574997`,14},{-4.0725`,1.3575`,6.7875`,14},{-4.0725`,6.7875`,-4.0725`,14},{-4.0725`,6.7875`,1.3575`,14},{-4.0725`,6.7875`,6.7875`,14},{1.3574997`,-4.0725`,-4.0725`,14},{1.3574997`,-4.0725`,1.3574997`,14},{1.3574997`,-4.0725`,6.7875`,14},{1.3575`,1.3575`,-4.0725`,14},{1.3575`,1.3575`,1.3575`,14},{1.3575`,1.3575`,6.7875`,14},{1.3575`,6.7875`,-4.0725`,14},{1.3575`,6.7875`,1.3575`,14},{1.3575`,6.7875`,6.7875`,14},{6.7875`,-4.0725`,-4.0725`,14},{6.7875`,-4.0725`,1.3575`,14},{6.7875`,-4.0725`,6.7875`,14},{6.7875`,1.3575`,-4.0725`,14},{6.7875`,1.3575`,1.3575`,14},{6.7875`,1.3575`,6.7875`,14},{6.7875`,6.7875`,-4.0725`,14},{6.7875`,6.7875`,1.3575`,14},{6.7875`,6.7875`,6.7875`,14},{-6.7875`,-1.3575`,1.3574995`,14},{-6.7875`,4.0725`,1.3574997`,14},{-1.3575`,-6.7875`,1.3574995`,14},{-1.3575`,-1.3575`,-4.0725`,14},{-1.3575`,-1.3575`,1.3574997`,14},{-1.3575`,-1.3575`,6.7875`,14},{-1.3574997`,4.0725`,-4.0725`,14},{-1.3574997`,4.0725`,1.3575`,14},{-1.3574997`,4.0725`,6.7875`,14},{-1.3574993`,9.5025`,1.3575`,14},{-1.3574993`,9.5025`,6.7875`,14},{4.0725`,-6.7875`,1.3574997`,14},{4.0725`,-1.3575`,-4.0725`,14},{4.0725`,-1.3575`,1.3575`,14},{4.0725`,-1.3575`,6.7875`,14},{4.0725`,4.0725`,-4.0725`,14},{4.0725`,4.0725`,1.3575`,14},{4.0725`,4.0725`,6.7875`,14},{4.0725`,9.5025`,-4.072499`,14},{4.0725`,9.5025`,1.3575`,14},{4.0725`,9.5025`,6.7875`,14},{9.5025`,-1.3575`,1.3575`,14},{9.5025`,-1.3575`,6.7875`,14},{9.5025`,4.0725`,-4.072499`,14},{9.5025`,4.0725`,1.3575`,14},{9.5025`,4.0725`,6.7875`,14},{-4.0725`,-1.3575`,-1.3575`,14},{-4.0725`,-1.3575`,4.0725`,14},{-4.0725`,4.0725`,-1.3575`,14},{-4.0725`,4.0725`,4.0725`,14},{-4.0725`,4.0725`,9.5025`,14},{-4.072499`,9.5025`,4.0725`,14},{1.3574995`,-6.7875`,-1.3575`,14},{1.3574995`,-6.7875`,4.0725`,14},{1.3575`,-1.3575`,-6.7875`,14},{1.3575`,-1.3575`,-1.3575`,14},{1.3575`,-1.3575`,4.0725`,14},{1.3575`,-1.3575`,9.5025`,14},{1.3575`,4.0725`,-6.7875`,14},{1.3575`,4.0725`,-1.3574997`,14},{1.3575`,4.0725`,4.0725`,14},{1.3575`,4.0725`,9.5025`,14},{1.3575`,9.5025`,-1.3574993`,14},{1.3575`,9.5025`,4.0725`,14},{6.7875`,-1.3575`,-1.3574995`,14},{6.7875`,-1.3575`,4.0725`,14},{6.7875`,-1.3575`,9.5025`,14},{6.7875`,4.0725`,-1.3574993`,14},{6.7875`,4.0725`,4.0725`,14},{6.7875`,4.0725`,9.5025`,14},{6.7875`,9.5025`,-1.3574988`,14},{6.7875`,9.5025`,4.0725`,14},{-6.7875`,1.3575`,-1.3575`,14},{-6.7875`,1.3575`,4.0725`,14},{-1.3575`,-4.0725`,-1.3575`,14},{-1.3575`,-4.0725`,4.0725`,14},{-1.3575`,1.3575`,-6.7875`,14},{-1.3575`,1.3575`,-1.3575`,14},{-1.3575`,1.3575`,4.0725`,14},{-1.3575`,1.3575`,9.5025`,14},{-1.3574995`,6.7875`,-1.3574995`,14},{-1.3574995`,6.7875`,4.0725`,14},{-1.3574995`,6.7875`,9.5025`,14},{4.0725`,-4.0725`,-1.3575`,14},{4.0725`,-4.0725`,4.0725`,14},{4.0725`,-4.0725`,9.5025`,14},{4.0725`,1.3575`,-6.7875`,14},{4.0725`,1.3575`,-1.3574997`,14},{4.0725`,1.3575`,4.0725`,14},{4.0725`,1.3575`,9.5025`,14},{4.0725`,6.7875`,-1.3574993`,14},{4.0725`,6.7875`,4.0725`,14},{4.0725`,6.7875`,9.5025`,14},{9.5025`,-4.0725`,4.0725`,14},{9.5025`,1.3575`,-1.3574993`,14},{9.5025`,1.3575`,4.0725`,14},{9.5025`,6.7875`,-1.3574988`,14},{9.5025`,6.7875`,4.0725`,14}};

2*Abs[Min[dccxyz4]] (* 14 Angstrom Diameter Silicon Sphere *)

Phase Contrast

α=.0005;mppx=((2.54 (10^(-2)))/(2400 240000));Δz=80*(10^(-10));Cs=1.2*(10^(-3));λ=0.0197*(10^(-10));Δf=-561*(10^(-10));

fieldWidth=N[npix*mppx] (* in meters *)

fieldWidth=20*(10^(-10)); (* in meters *)

mag=((2.54*(10^(-2)))/(2400*fieldWidth/npix)) (* magnification of image digitized at 2400 dpi *)

metersPerPixel[fw_]:=(fw/npix)

AngstromsPerPixel[fwA_]:=(fwA/npix)

freqNyquist=N[(npix/(2*fieldWidth))] (* in reciprocal Angstroms *)

1+1

IncidentBeam=Table[1.0,{x,1,npix},{y,1,npix}];

exitSurfacePhase[projectedPotential_]:=(* MapThread[Times,{FLA, *)Map[Exp,-I projectedPotential,{2}](*,IncidentBeam},2] *);

CircleAp[x_,y_,r_]:=If[Sqrt[(x^2)+(y^2)]<r,1.0,0.1]

OA=aOff[Table[CircleAp[u-nh,v-nh,nh],{u,1,npix},{v,1,npix}]];

backFocalPlane[exitSurfacePhase_]:=Fourier[N[exitSurfacePhase]];

objectiveApertureWave[backFocalPlane_]:=MapThread[Times,{OA,backFocalPlane},2];

ω[u_,v_]:=(Sqrt[(((((u-nh)/fieldWidth)))^2)+(((((v-nh)/fieldWidth)))^2)]);

χ[ω_]:=π Δf λ (ω^2)+(1/2)π Cs (λ^3)(ω^4)

Phase2D[A_,{u_,v_}]:=A (E^(I χ[ω[u,v]]))

PhaseChange=MapIndexed[Phase2D,Table[1,{npix},{npix}],{2}];

Damp2D[Δz_,{u_,v_}]:=(E^(-(((π λ Δz))^2) ((ω[u,v])^4)/2))

Damp=MapIndexed[Damp2D,Table[Δz,{npix},{npix}],{2}];

Conv2D[α_,{u_,v_}]:=(E^(-(((((π α)/λ)))^2) (((Cs (λ^3) ((ω[u,v])^3)+Δf λ ω[u,v]))^2)))

ConvMatrix=MapIndexed[Conv2D,Table[α,{npix},{npix}],{2}];

Damp=MapThread[Times,{Damp,ConvMatrix},2];

CTF=aOff[MapThread[Times,{Damp,PhaseChange},2]]; (* presently CTF does not depend on fieldWidth in Meters which is fixed *)

ωnew[u_,v_,fwm_]:=(Sqrt[(((((u-nh)/fwm)))^2)+(((((v-nh)/fwm)))^2)])

CTFnew[fwm_,Δf_]:=aOff[Table[(E^(-(((((π α)/λ)))^2) (((Cs (λ^3) ((ωnew[u,v,fwm])^3)+Δf λ ωnew[u,v,fwm]))^2)))(E^(-(((π λ Δz))^2) (((ωnew[u,v,fwm])^4)/2)))(E^(I χ[ωnew[u,v,fwm]])),{u,1,npix},{v,1,npix}]]

fieldWidth

CTFnew[fieldWidth,Δf][[1,1]]

Dimensions[CTF]

CTF[[1,1]]

imageAmplitude[objectiveApertureWave_]:=InverseFourier[MapThread[Times,{CTF,objectiveApertureWave},2]];

ListDensityPlot[ArrayPwr[imageAmplitude[]],Mesh->False,PlotRange->{-1,2},Frame->False,ColorFunction->ColorData["GrayTones"],ImageSize->npix]

Utilities

zonexyz[u_,v_,w_]:={{5.43,0,0},{0,5.43,0},{0,0,5.43}}.{u,v,w} (* this gives the viewing direction in cartesian coordinates? *)

tilt[zone_]:=If[(zone[[1]]==0)&&(zone[[2]]==0),0,ArcTan[zone[[3]],Norm[{zone[[1]],zone[[2]]}]]] (* this determines the corresponding tilt *)

azimuth[zone_]:=If[(zone[[1]]==0)&&(zone[[2]]==0),0,ArcTan[zone[[1]],zone[[2]]]] (* this determines the corresponding azimuth *)

zonexyz[0,0,1]

{tilt[zonexyz[0,0,1]],azimuth[zonexyz[0,0,1]]}*180/Pi

RotationMatrix[tilt[zonexyz[1,1,1]],{Cos[azimuth[zonexyz[1,1,1]]],Sin[azimuth[zonexyz[1,1,1]]],0}]//MatrixForm

RotationMatrix[tilt[zonexyz[0,0,1]],{Cos[azimuth[zonexyz[0,0,1]]],Sin[azimuth[zonexyz[0,0,1]]],0}]//MatrixForm

RotationMatrix[-ϕ,{0,0,1}]

Range[5]

ProjPot4x[npixels_, molecule_, fieldwidth_,rotate_,tilt_,fieldoffset_]:=
(* Here molecule is a list of atom positions of the form {{x1,y1,z1,Z1},...,{xN,yN,zN,ZN}}, 
npix are the number of pixels on a side in the output array, fieldwidth is the width in x and y that the output array spans, beameuler is so far not used, and fieldoffset relative to the origin of molecule coordinates has the form {xoff,yoff} *)
Module[
{scratch=Table[0,{i,1,npixels},{j,1,npixels}],mrot=Table[{0,0,0},{i,1,Length[molecule]}],indx,nndx,Fxy,infield,kk},
(* begin local declarations *)
For[i=1,i<Length[molecule],i++,
mrot[[i]]= RotationMatrix[tilt,{1,0,0}].RotationMatrix[rotate,{0,0,1}].({molecule[[i,1]],molecule[[i,2]],molecule[[i,3]]}-fieldoffset-(fieldwidth/2){1,1,0})+(fieldwidth/2){1,1,0} (* this gets the rotated coordinate set **)
];(*  Print[mrot]; *)
(* scratch=Table[0,{i,1,npixels},{j,1,npixels}];  *)
(* this creates a scratch table filled with zeros *) 
indx[m_,n_]:=Round[((mrot[[m,n]])/fieldwidth)*npixels+1];
 (* calculate the main indx[m,n]: the n=1,2 pixel index associated with the mth atom in molecule. *)
nndx[m_,n_]:=
If[((mrot[[m,n]])/fieldwidth)*npixels+1-indx[m,n]>0,indx[m,n]+1,indx[m,n]-1];
(* nndx[m,n]: the diagonal nearest-neighbor index in the nth direction associated with the mth atom *)
Fxy[m_,n_]:=1-Abs[((mrot[[m,n]])/fieldwidth)*npixels+1-indx[m,n]];
(* fraction Fxy[m,n] is Fx if n=1 and Fy if n=2 frac btwn {1/2,1} for the mth atom in molecule. *)
infield[m_]:=If[(indx[m,1]≥2)&&(indx[m,1]≤npixels-1)&&(indx[m,2]≥2)&&(indx[m,2]≤npixels-1),True,False]; (* Now apply: record atom location in scratch if it's in the field of view *)
Do[ (* now cycle through all atoms *)
	If[infield[k], (* if and only if the kth atom is in the field of view... *)
	Part[scratch,indx[k,1],indx[k,2]]=Part[scratch,indx[k,1],indx[k,2]]+molecule[[k,4]] Fxy[k,1] Fxy[k,2]];

	If[infield[k],Part[scratch,nndx[k,1],indx[k,2]]=
Part[scratch,nndx[k,1],indx[k,2]]+molecule[[k,4]] (1-Fxy[k,1]) Fxy[k,2]];

	If[infield[k],Part[scratch,indx[k,1],nndx[k,2]]=
Part[scratch,indx[k,1],nndx[k,2]]+molecule[[k,4]] Fxy[k,1] (1-Fxy[k,2])];

	If[infield[k],Part[scratch,nndx[k,1],nndx[k,2]]=
Part[scratch,nndx[k,1],nndx[k,2]]+molecule[[k,4]] (1-Fxy[k,1]) (1-Fxy[k,2])],
	(* first add lion's share of atom type [[k,4]] to scratch array at location indx[k,1],indx[k,2]; the next 3 lines distribute "spillover" using nndx and Fxy to increment 3 neighboring cells *)
	{k,1,Dimensions[molecule,1][[1]]}
];
N[ scratch]  (* finally, return the scratch array value, where N[] has been added to gloss over an error? *) 
]

viewhalfangle[fw_]:=fw*Pi/30000; viewvectormag=10000;

N[viewhalfangle[19.5]]

N[14/10000]

N[((14/2)/10000)]

viewHalfAngle2[fw_]:=(21/7)ArcTan[10000,(fw/2)] (* factor of 21/7 is of mysterious origin *)

N[viewHalfAngle2[14]]

Degree

35 Degree

N[Sin[35 Degree]]

N[Sin[35*Pi/180]]

Interactive Module

Manipulate[
GraphicsGrid[
{
{
Graphics3D[{DotsZ[dccxyz4,0,0,0,.5*25/fw](*, Sphere[{0,0,0},5] *)},Boxed->False, 
ViewVector ->  {
10000* RotationMatrix[(*tilt amount:*)θview*(π/180),(*tilt axis:*){Cos[ϕview*(π/180)],Sin[ϕview*(π/180)],0}].{0,0,1}+{-offset[[2]]+2.5,-offset[[1]]+2.5,2.5},
 {offset[[1]]+2.5,offset[[2]]+2.5,2.5} }
(*
ViewPoint->600* RotationMatrix[tilt[zonexyz[u,v,w]],{Cos[azimuth[zonexyz[u,v,w]]],Sin[azimuth[zonexyz[u,v,w]]],0}].{0,0,1}(* 100*zonexyz[u,v,w]*) *),ViewVertical-> 600*RotationMatrix[(*tilt amount:*)θview*(π/180),(*tilt axis:*){Cos[ϕview*(π/180)],Sin[ϕview*(π/180)],0}].(*vertical direction before tilt:*){-Sin[ϕview*(π/180)],Cos[ϕview*(π/180)],0},
ViewAngle->fw*Pi/30000,AspectRatio->1, ImageSize->128 , RotationAction->"Clip"] ,
If[mode==1,
ListDensityPlot[
Transpose[ProjPot4x[npix,dccxyz4,fw,ϕview*(π/180),θview*(π/180),{offset[[1]]-fw/2+2.5,offset[[2]]-fw/2+2.5,2.5}]],
Frame->False ,ColorFunction->(Blend[Fire5,#]&),ImageSize->128
],
If[mode==2,
ListDensityPlot[
aOff[
Abs[
Fourier[
Transpose[
ProjPot4x[npix,dccxyz4,fw,ϕview*(π/180),θview*(π/180),{offset[[1]]-fw/2+2.5,offset[[2]]-fw/2+2.5,2.5}]
]
]
]
],
Frame->False,ColorFunction->(Blend[Fire5,#]&), ImageSize->128
],
If[mode==3,
CTF=CTFnew[fw*(10^(-10)),Δf];
ListDensityPlot[
ArrayPwr[
imageAmplitude[
objectiveApertureWave[
backFocalPlane[
exitSurfacePhase[
Transpose[
ProjPot4x[npix,dccxyz4,fw,ϕview*(π/180),θview*(π/180),{offset[[1]]-fw/2+2.5,offset[[2]]-fw/2+2.5,2.5}]
]
]
]
]
]
],
Mesh->False,
PlotRange->{-1,2},
Frame->False,
ColorFunction->ColorData["GrayTones"],
ImageSize->npix
],
Graphics[{White,Rectangle[]},ImageSize->128]
]
]
]
}
}
],
{{fw,50,"field width in Å"},1,50,Appearance->"Labeled"},
{{offset,{0,0},"offset in Å"},{-15,-15},{15,15},Appearance->"Labeled"},
(* {{direction,{0,0},"rotate-tilt"},{-Pi,-Pi/2},{Pi,Pi/2},Appearance->"Labeled"}, *)
{{θview,0,"tilt"},0,90,Appearance->"Labeled"},
{{ϕview,0,"rotate"},-180,180,Appearance->"Labeled"},
{{mode,0,"mode"},{0->"NA",1-> "PrjPot",2->"DifPat",3->"BriFld",4->"PowSpc",5->"DrkFld"},ControlType->SetterBar},
(* NA means not applicable, PP means projected potential, DP means diffraction pattern, BF means brightfield image, PS means power spectrum of image, DF means darkfield image, ... *)
(* 
{{u,0},-4,4,Appearance->"Labeled"},
{{v,0},-4,4,Appearance->"Labeled"},
{{w,1},-4,4,Appearance->"Labeled"},
*)
(* "test", *)
Deployed->True,
SaveDefinitions->True,
ContinuousAction->False
]