'------------------------------------------------------------------------------------------------------------------------------------- 'Useful macros collected from several sources not included in TOPAS.INC '------------------------------------------------------------------------------------------------------------------------------------- 'Stephens model of anisotropic peak broadening 'The following are taken from http://topas.dur.ac.uk/topaswiki/doku.php?id=stephens_peakshape&s[]=stephens and relate to the Stephens model of anistropic peak broadening 'P.W. Stephens, J. Appl. Cryst. (1999) 32, 281-9 'Stephens_tetragonal(s400, s004, s220, s202, eta) macro Stephens_tetragonal(s400, s004, s220, s202, eta) { local mhkl = Abs(s400 (H^4+ K^4)+ s004 L^4+ s220 H^2 K^2 + s202 (H^2 L^2 + K^2 L^2) ); local pp = D_spacing^2 * Sqrt(Max(mhkl,0)) / 1000; gauss_fwhm = 1.8/3.1415927 pp (1-eta) Tan(Th) + 0.0001; lor_fwhm = 1.8/3.1415927 pp eta Tan(Th) + 0.0001; } 'Stephens_monoclinic(s400, s040, s004, s220, s202, s022, s301, s121, s103, eta) macro Stephens_monoclinic(s400, s040, s004, s220, s202, s022, s301, s121, s103, eta) { local mhkl = H^4 s400 + K^4 s040 + L^4 s004 + H^2 K^2 s220 + H^2 L^2 s202 + K^2 L^2 s022 + H K^2 L s121 + H L^3 s103 + H^3 L s301; local pp = D_spacing^2 * Sqrt(Max(mhkl,0)) / 1000; gauss_fwhm = 1.8/3.1415927 pp (1-eta) Tan(Th) + 0.0001; lor_fwhm = 1.8/3.1415927 pp eta Tan(Th) + 0.0001; } 'Stephens_hexagonal(s400, s202, s004, eta) macro Stephens_hexagonal(s400, s202, s004, eta) { local mhkl = H^4 s400 + K^4 s400 + L^4 s004 + H^2 K^2 3 s400 + H^2 L^2 s202 + K^2 L^2 s202 + H K L^2 s202 + H^3 K 2 s400 + H K^3 2 s400; local pp = D_spacing^2 * Sqrt(Max(mhkl,0)) / 1000; gauss_fwhm = 1.8/3.1415927 pp (1-eta) Tan(Th) + 0.0001; lor_fwhm = 1.8/3.1415927 pp eta Tan(Th) + 0.0001; } 'Stephens_orthorhombic(s400, s040, s004, s220, s202, s022, eta) macro Stephens_orthorhombic(s400, s040, s004, s220, s202, s022, eta) { local mhkl = H^4 s400 + K^4 s040 + L^4 s004 + H^2 K^2 s220 + H^2 L^2 s202 + K^2 L^2 s022; local pp = D_spacing^2 * Sqrt(Max(mhkl,0)) / 1000; gauss_fwhm = 1.8/3.1415927 pp (1-eta) Tan(Th) + 0.0001; lor_fwhm = 1.8/3.1415927 pp eta Tan(Th) + 0.0001; } 'To calculate a pattern macro calc_pattern(start,stop,step,filename) {yobs_eqn !aac.xy = 1; min start max stop del step xdd_out filename load out_record out_fmt out_eqn { " %11.5f " = X; " %11.5f\n" = Ycalc; } } 'To output a Reitveld plot macro Rietveld_Plot(filename) { xdd_out filename load out_record out_fmt out_eqn { " %11.5f " = X; " %11.5f " = Yobs; " %11.5f " = Ycalc; " %11.5f\n" = Yobs-Ycalc; } } 'Outputs a cif for use with crystalmaker 'Taken from http://topas.dur.ac.uk/topaswiki/doku.php?id=out_cif_crystalmaker macro Out_CIF_crystalmaker(file) {Out_CIF_crystalmaker_1(file) Out_CIF_crystalmaker_2(file)} macro Out_CIF_crystalmaker_1(file) { out file Out_String("\ndata_") Out(Get(phase_name), "\n_chemical_name_mineral ?%s?") Out(Get(a), "\n_cell_length_a %V") Out(Get(b), "\n_cell_length_b %V") Out(Get(c), "\n_cell_length_c %V") Out(Get(al), "\n_cell_angle_alpha %V") Out(Get(be), "\n_cell_angle_beta %V") Out(Get(ga), "\n_cell_angle_gamma %V") Out(Get(cell_volume), "\n_cell_volume %V") Out(Get(sp_grp_char), "\n_space_group %s") Out_String("\nloop_\n_symmetry_equiv_pos_as_xyz") Out(Get(sp_xyzs_txt), "%s") Out_String("\nloop_") Out_String("\n_atom_site_label") ' Out_String("\n_atom_site_type_symbol") ' Out_String("\n_atom_site_symmetry_multiplicity") Out_String("\n_atom_site_fract_x") Out_String("\n_atom_site_fract_y") Out_String("\n_atom_site_fract_z") Out_String("\n_atom_site_occupancy") Out_String("\n_atom_site_B_iso_or_equiv") atom_out file append load out_record out_fmt out_eqn { "\n%s" = Get_From_String(Get(current_atom), site); ' " %s" = Get_From_String(Get(current_atom), atom); ' " %3.0f" = Get_From_String(Get(current_atom), num_posns); " %V" = Get_From_String(Get(current_atom), x); " %V" = Get_From_String(Get(current_atom), y); " %V" = Get_From_String(Get(current_atom), z); " %V" = Get_From_String(Get(current_atom), occ); " %V" = Get_From_String(Get(current_atom), beq); } } macro Out_CIF_crystalmaker_2(file) { out file append Out_String("\nloop_") Out_String("\n_atom_site_aniso_label") ' Out_String("\n_atom_site_type_symbol") Out_String("\n_atom_site_aniso_U_11") Out_String("\n_atom_site_aniso_U_22") Out_String("\n_atom_site_aniso_U_33") Out_String("\n_atom_site_aniso_U_12") Out_String("\n_atom_site_aniso_U_13") Out_String("\n_atom_site_aniso_U_23") atom_out file append load out_record out_fmt out_eqn { "\n%s" = Get_From_String(Get(current_atom), site); ' " %s" = Get_From_String(Get(current_atom), atom); " %V" = Get_From_String(Get(current_atom), u11); " %V" = Get_From_String(Get(current_atom), u22); " %V" = Get_From_String(Get(current_atom), u33); " %V" = Get_From_String(Get(current_atom), u12); " %V" = Get_From_String(Get(current_atom), u13); " %V" = Get_From_String(Get(current_atom), u23); } } macro Fix_Biso(p2,v2) { prm p2 v2 min -10 max 10 scale_pks = (Exp(-2*p2 Sin(Th)^2/Lam^2)) ; } macro wifd_pkshape(two_theta,two_thetav,aprm,aprmv,tauf_2,tauf_2v,tauf_1,tauf_1v,tauf_0,tauf_0v,taus_0,taus_0v,taus_1,taus_1v,t_eff,t_effv,hhh,hhhv,hhh_0,hhh_0v) { local !two_theta two_thetav local !lam = 2 D_spacing Sin(0.00872664626 two_theta) ; prm aprm aprmv min 1.5 max 2.5 prm tauf_2 tauf_2v min 0.01 prm tauf_1 tauf_1v min 0.01 prm tauf_0 tauf_0v min 0.001 user_defined_convolution = X^aprm Exp(-X / (tauf_0 + tauf_1 lam + tauf_2 lam^2)); min 0 max = 20 aprm tauf_1 lam; prm taus_0 taus_0v min 0.001 prm taus_1 taus_1v min 0.001 prm t_eff t_effv min 90.0 local storage = Exp( -950 /(t_eff lam lam)) ; push_peak prm hhh hhhv min .001 prm hhh_0 hhh_0v min 0.001 hat = hhh_0 + hhh D_spacing; scale_top_peak = (1-storage); bring_2nd_peak_to_top hat = hhh_0 + hhh D_spacing; scale_top_peak = storage; exp_conv_const = -Ln(0.001) (taus_0 + lam taus_1); add_pop_1st_2nd_peak } macro wifd_mic(two_theta,two_thetav,tspow,tspowv,tauf_0,tauf_0v,tauf_1,tauf_1v,taus_0,taus_0v,taus_1,taus_1v,t_eff,t_effv,hhh,hhhv,double_pulse,double_pulsev) { local !two_theta two_thetav prm tspow tspowv min 1.5 max 5.5 prm tauf_0 tauf_0v min 0.001 prm tauf_1 tauf_1v min 0.01 prm taus_0 taus_0v min 0.001 prm taus_1 taus_1v min 0.001 prm t_eff t_effv min 90.0 prm hhh hhhv min .001 prm double_pulse double_pulsev min 0.001 local !lam = 2 D_spacing Sin(0.00872664626 two_theta) ; user_defined_convolution = X^tspow Exp(-X / (tauf_0 + tauf_1 lam )); min 0 max = 10 (tauf_0 + tauf_1 lam ); local !storage = Exp( -950 /(t_eff lam^2)) hhh; push_peak hat = double_pulse ; scale_top_peak = (1-storage)/double_pulse; bring_2nd_peak_to_top exp_conv_const = -Ln(0.001) / (taus_0 + taus_1 / lam^2); scale_top_peak = storage/double_pulse; add_pop_1st_2nd_peak } macro wifd_mic_new(two_theta,two_thetav,tspow,tspowv,tauf_0,tauf_0v,tauf_1,tauf_1v,taus_0,taus_0v,taus_1,taus_1v,t_eff,t_effv,hhh,hhhv,double_pulse,double_pulsev) { local !two_theta two_thetav prm tspow tspowv min 1.5 max 5.5 prm tauf_0 tauf_0v min 0.001 prm tauf_1 tauf_1v min 0.01 prm taus_0 taus_0v min 0.001 prm taus_1 taus_1v min 0.001 prm t_eff t_effv min 90.0 prm hhh hhhv min .001 prm double_pulse double_pulsev min 0.001 local !lam = 2 D_spacing Sin(0.00872664626 two_theta) ; user_defined_convolution = X^tspow Exp(-X / (tauf_0 + tauf_1 lam )); min 0 max = 10 (tauf_0 + tauf_1 lam ); local !storage = Exp( -950 /(t_eff lam^2)) hhh; push_peak hat = double_pulse ; scale_top_peak = (1-storage); bring_2nd_peak_to_top exp_conv_const = -Ln(0.001) / (taus_0 + taus_1 / lam^2); scale_top_peak = storage; add_pop_1st_2nd_peak } macro tof_sample_peakshape(lor,lor_val,dsp,dsp_val,dspsq,dspsq_val) { prm dsp dsp_val del = 0.05 Val + 1; min 1 prm dspsq dspsq_val del = 0.05 Val + 1; min 0 peak_type pv pv_lor lor lor_val pv_fwhm = dsp D_spacing + dspsq D_spacing^2; } macro synchrotron_peakshape(pr1, pr1val, pr2, pr2val, pr3, pr3val, pr4, pr4val) { prm pr1 pr1val min 0 prm pr2 pr2val min 0 prm pr3 pr3val min 0 prm pr4 pr4val min 0 lor_fwhm = pr1 Tan(Th) + pr2/Cos(Th) ; gauss_fwhm = pr3 Tan(Th) + pr4/Cos(Th) ; } macro inverse_Axial_Model(c, v) { #m_argu c If_Prm_Eqn_Rpt(c, v, min 0.0001 max 50) circles_conv = 0.5 Rad ( CeV(c,v) / Rs)^2 / Tan(2 Th); } macro output_tickmarks(filename) { phase_out filename load out_record out_fmt out_eqn { " %11.5f" = D_spacing; " %11.5f" = 2 Rad Th; " %3.0f\n" = 0; } } macro bkg_lognormal(scl, sclv, mu, muv, sig, sigv) { #m_argu scl #m_argu mu #m_argu sig prm scl sclv prm mu muv prm sig sigv fit_obj = scl Exp((-(Ln(X) - mu)^2)/(2*sig^2)) 1/(X*Sqrt(2*Pi)*sig); } macro bkg_Gumbel2(scl, sclv, para, parav, parb, parbv) { #m_argu scl #m_argu para #m_argu parb prm scl sclv prm para parav prm parb parbv fit_obj = scl para parb X^(-para-1) Exp(-parb X^(-para)); } macro bkg_hlogd(scl, sclv, parak, parakv) { #m_argu scl #m_argu parak prm scl sclv prm parak parakv fit_obj = scl (2 Exp(-parak X))/((1+Exp(-parak X))^2); } macro bkg_hill(Vmax, Vmaxv, Km, Kmv, paran, paranv) { #m_argu Vmax #m_argu Km #m_argu paran prm Vmax Vmaxv prm Km Kmv prm paran paranv fit_obj = (Vmax*X^paran)/(Km^paran + X^paran); } macro One_on_X_power(Vmax, Vmaxv, pow, powv, xstart, xstartv) { #m_argu Vmax #m_argu pow #m_argu xstart prm Vmax Vmaxv prm pow powv prm xstart xstartv fit_obj = vmax / (X-xstart)^pow; } macro AnisoCS( mod, h11, k11, l11, h22, k22, l22, ac, av, bc, bv, cc, cv, pc, pv, tc, tv, nc, nv, sp) { ' As published in Ectors et al. (2015) J. Appl. Cryst. 48, 189-194 ' mod which model to use 1 for ellipsoid 2 for elliptic cylinder 3 for cuboid ' h11 k11 l11 vector defining the z-axis ' h22 k22 l22 vector defining the x-axis ' av bv cv rx ry rz of the model ' pv tv nv additional rotation parameters theta' theta' nu ' sp activate special method for triclinic systems 1/0 on/off ' --------- Code starts here ------ peak_buffer_step 0 'Initialize variables #m_argu mod #m_argu ac If_Prm_Eqn_Rpt(ac, av, min 0.0001 max 10000) #m_argu bc If_Prm_Eqn_Rpt(bc, bv, min 0.0001 max 10000) #m_argu cc If_Prm_Eqn_Rpt(cc, cv, min 0.0001 max 10000) #m_argu pc If_Prm_Eqn_Rpt(pc, pv, min 0 max 180) #m_argu tc If_Prm_Eqn_Rpt(tc, tv, min 0 max 180) #m_argu nc If_Prm_Eqn_Rpt(nc, nv, min 0 max 360) 'G-Metrictensor just for output local !m11 = Get(a)^2; local !m12 = Get(a)*Get(b)*Cos(Get(ga)*Deg); local !m13 = Get(a)*Get(c)*Cos(Get(be)*Deg); local !m21 = m12; local !m22 = Get(b)^2; local !m23 = Get(b)*Get(c)*Cos(Get(al)*Deg); local !m31 = m13; local !m32 = m23; local !m33 = Get(c)^2; 'G*-Reciprocal metrictensor and reciprocal lattice parameters local !vol = Get(a)*Get(b)*Get(c)*Sqrt(1-((Cos(Get(al)*Deg))^2)-((Cos(Get(be)*Deg))^2)-((Cos(Get(ga)*Deg))^2)+2*(Cos(Get(al)*Deg)*Cos(Get(be)*Deg)*Cos(Get(ga)*Deg))); local !ar = (Get(b)*Get(c)*Sin(Get(al)*Deg))/vol; local !br = (Get(a)*Get(c)*Sin(Get(be)*Deg))/vol; local !cr = (Get(b)*Get(a)*Sin(Get(ga)*Deg))/vol; local !alr = (Cos(Get(be)*Deg)*Cos(Get(ga)*Deg)-Cos(Get(al)*Deg))/(Sin(Get(be)*Deg)*Sin(Get(ga)*Deg)); local !ber = (Cos(Get(al)*Deg)*Cos(Get(ga)*Deg)-Cos(Get(be)*Deg))/(Sin(Get(al)*Deg)*Sin(Get(ga)*Deg)); local !gar = (Cos(Get(be)*Deg)*Cos(Get(al)*Deg)-Cos(Get(ga)*Deg))/(Sin(Get(be)*Deg)*Sin(Get(al)*Deg)); local !r11 = ar^2; local !r12 = ar*br*gar; local !r13 = ar*cr*ber; local !r21 = r12; local !r22 = br^2; local !r23 = br*cr*alr; local !r31 = r13; local !r32 = r23; local !r33 = cr^2; 'Test if special method triclinic is activated #m_ifarg sp 1 local !h1 = 0; local !h2 = If(alr <= 0.00001, 0 , 1); local !k1 = 0; local !k2 = If(alr <= 0.00001, 1 , -(r31 / r32)); local !l1 = 1; local !l2 = 0; #m_else local !h1 = h11; local !h2 = h22; local !k1 = k11; local !k2 = k22; local !l1 = l11; local !l2 = l22; #m_endif '|hkl|Prim aka z-axis local rowh1 = (h1 r11 + k1 r12 + l1 r13)*h1; local rowk1 = (h1 r21 + k1 r22 + l1 r23)*k1; local rowl1 = (h1 r31 + k1 r32 + l1 r33)*l1; local dvec1 = Sqrt(rowh1+rowk1+rowl1); '|hkl|Sec aka x-axis local rowh2 = (h2 r11 + k2 r12 + l2 r13)*h2; local rowk2 = (h2 r21 + k2 r22 + l2 r23)*k2; local rowl2 = (h2 r31 + k2 r32 + l2 r33)*l2; local dvec2 = Sqrt(rowh2+rowk2+rowl2); '|hkl|Current local rowH = (H r11 + K r12 + L r13)*H; local rowK = (H r21 + K r22 + L r23)*K; local rowL = (H r31 + K r32 + L r33)*L; local dvec3 = Sqrt(rowH+rowK+rowL); 'Angle phi local rowg1 = (H r11 + K r12 + L r13)*h1; local rowg2 = (H r21 + K r22 + L r23)*k1; local rowg3 = (H r31 + K r32 + L r33)*l1; local sum1 = rowg1+rowg2+rowg3; local cosphi = (sum1/(dvec1*dvec3)); local sinphi = Sqrt(Abs(1-(cosphi^2))); 'Angle theta local rowp1 = (H r11 + K r12 + L r13)*h2; local rowp2 = (H r21 + K r22 + L r23)*k2; local rowp3 = (H r31 + K r32 + L r33)*l2; local sum2 = rowp1+rowp2+rowp3; local coseta =(sum2/(dvec2*dvec3)); local costheta = If(sinphi <= 0.00001, 1, coseta/sinphi); local sintheta = Sqrt(Abs(1-(costheta^2))); 'Handedness local p = h1*(k2*L - l2*K) + k1*(l2*H - h2*L) + l1*(h2*K - k2*H); local pre = If (p <= 0 , -1 , 1); 'Rotation: Rotationmatrix and calculation of phirot and thetarot local grot11 = (Cos(CeV(tc,tv) * Deg)*Cos(CeV(pc,pv) * Deg)*Cos(CeV(nc,nv) * Deg))-(Sin(CeV(tc,tv) * Deg)*Sin(CeV(nc,nv) * Deg)); local grot21 = -(Cos(CeV(tc,tv) * Deg)*Cos(CeV(pc,pv) * Deg)*Sin(CeV(nc,nv) * Deg))-(Sin(CeV(tc,tv) * Deg)*Cos(CeV(nc,nv) * Deg)); local grot31 = Cos(CeV(tc,tv) * Deg)*Sin(CeV(pc,pv) * Deg); local grot12 = (Sin(CeV(tc,tv) * Deg)*Cos(CeV(pc,pv) * Deg)*Cos(CeV(nc,nv) * Deg))+(Cos(CeV(tc,tv) * Deg)*Sin(CeV(nc,nv) * Deg)); local grot22 = -(Sin(CeV(tc,tv) * Deg)*Cos(CeV(pc,pv) * Deg)*Sin(CeV(nc,nv) * Deg))+(Cos(CeV(tc,tv) * Deg)*Cos(CeV(nc,nv) * Deg)); local grot32 = Sin(CeV(tc,tv) * Deg)*Sin(CeV(pc,pv) * Deg); local grot13 = -(Sin(CeV(pc,pv) * Deg)*Cos(CeV(nc,nv) * Deg)); local grot23 = Sin(CeV(pc,pv) * Deg)*Sin(CeV(nc,nv) * Deg); local grot33 = Cos(CeV(pc,pv) * Deg); local cosphirot = (costheta*sinphi*grot31)+(pre*sintheta*sinphi*grot32)+(cosphi*grot33); local sinphirot = Sqrt(Abs(1-(cosphirot^2))); local xfac = (costheta*sinphi*grot11)+(pre*sintheta*sinphi*grot12)+(cosphi*grot13); local costhetarot = If(sinphirot <= 0.0001, 1, xfac / sinphirot); local sinthetarot = Sqrt(Abs(1-(costhetarot^2))); 'Models: 'Ellipsoid #m_ifarg mod 1 local volop = ((4/3) 3.14159265358979 CeV(ac,av) CeV(bc,bv) CeV(cc,cv)); local areaop = 3.14159265358979 Sqrt(((CeV(cc,cv) sinphirot)^2)(((CeV(ac,av) sinthetarot)^2)+(CeV(bc,bv) costhetarot)^2)+((CeV(ac,av) CeV(bc,bv) cosphirot)^2)); local sizeL = volop / areaop; #m_endif 'Cylinder #m_ifarg mod 2 local volop = (3.14159265358979 CeV(ac,av) CeV(bc,bv) (2 * CeV(cc,cv))); local areaop = (Sqrt( ((CeV(ac,av) sinthetarot)^2) + ((CeV(bc,bv) costhetarot)^2)) 4 * CeV(cc,cv) * Abs(sinphirot)) + (3.14159265358979 CeV(ac,av) CeV(bc,bv) Abs(cosphirot)); local sizeL = volop / areaop; #m_endif 'Cuboid #m_ifarg mod 3 local volop = (8 CeV(ac,av) CeV(bc,bv) CeV(cc,cv)); local areaop = 4 CeV(ac,av) CeV(bc,bv) Abs(cosphirot) + 4 CeV(ac,av) CeV(cc,cv) Abs(sinphirot) Abs(sinthetarot) + 4 CeV(bc,bv) CeV(cc,cv) Abs(costhetarot) Abs(sinphirot) ; local sizeL = volop / areaop; #m_endif 'For Output purposes only local ax = CeV(ac,av); local bx = CeV(bc,bv); local cx = CeV(cc,cv); local px = CeV(pc,pv); local tx = CeV(tc,tv); local nx = CeV(nc,nv); local modx = mod; local vole = volop; local larea = volop^(1/3); local kk = larea / sizeL ; local areaoe = areaop; local phi = If(cosphi < -1, 180, If (cosphi > 1, 0, ArcCos(cosphi)*Rad)); local the = If(costheta < -1, 180, If (costheta > 1, 0, ArcCos(costheta)*Rad)); local phirot = If(cosphirot < -1, 180, If (cosphirot > 1, 0, ArcCos(cosphirot)*Rad)); local therot = If(costhetarot < -1, 180, If (costhetarot > 1, 0, ArcCos(costhetarot)*Rad)); 'CS_Macro 'lor_fwhm = 0.1 180/pi^2 * Lambda / Cos(Th) * sizeL ... 0.1 for nm->Angstroem lor_fwhm = 0.1 18.2378130556208 Lam / (Cos(Th) sizeL); } macro AnisoCSout(file) { ' As published in Ectors et al. (2015) J. Appl. Cryst. 48, 189-194 out file Out_String("Analysis Report: ") Out(Get(phase_name),"%s \n") Out(Get (r_wp), " Rwp: %8.3f \n") Out(Get (r_exp), " Rexp:%8.3f \n\n") Out_String("G-Metrictensor \n\n") Out(m11,"%8.3f") Out(m12,"%8.3f") Out(m13,"%8.3f\n") Out(m21,"%8.3f") Out(m22,"%8.3f") Out(m23,"%8.3f\n") Out(m31,"%8.3f") Out(m32,"%8.3f") Out(m33,"%8.3f\n\n") Out_String("G*-Metrictensor \n\n") Out(r11,"%8.3f") Out(r12,"%8.3f") Out(r13,"%8.3f\n") Out(r21,"%8.3f") Out(r22,"%8.3f") Out(r23,"%8.3f\n") Out(r31,"%8.3f") Out(r32,"%8.3f") Out(r33,"%8.3f\n\n") Out_String("Rotationmatrix \n\n") Out(grot11,"%8.3f") Out(grot12,"%8.3f") Out(grot13,"%8.3f\n") Out(grot21,"%8.3f") Out(grot22,"%8.3f") Out(grot23,"%8.3f\n") Out(grot31,"%8.3f") Out(grot32,"%8.3f") Out(grot33,"%8.3f\n\n") Out_String("Model parameters \n") Out_String(" (1) Ellipsoid / (2) Cylinder / (3) Cuboid \n") Out(modx, " Model: %1.0f \n\n") Out(h1," z-axis hkl: %8.3f") Out(k1," %8.3f") Out(l1," %8.3f \n") Out(px, " Rotation phi: %8.3f") Out(tx, " theta: %8.3f \n\n") Out(h2," x-axis hkl: %8.3f") Out(k2," %8.3f") Out(l2," %8.3f \n") Out(nx, " Additional rotation nu: %8.3f \n\n") Out(ax," rx-radius: %8.3f nm \n") Out(bx," ry-radius: %8.3f nm \n") Out(cx," rz-radius: %8.3f nm \n") Out(vole," Volume: %8.3f nm^3 \n") Out(larea," : %8.3f nm \n\n") Out_String(" H K L M 2Th D phi theta phirot thetarot \n") phase_out file append load out_record out_fmt out_eqn { "%4.0f" = H; "%4.0f" = K; "%4.0f" = L; "%4.0f" = M; "%8.3f" = 2 Th Rad; "%8.3f" = D_spacing; "%8.3f" = phi; "%8.3f" = the; "%8.3f" = phirot; "%8.3f" = therot; "%8.3f\n" = sizeL; } } macro AnisoCSg( tauc, tauv) { #m_argu tauc If_Prm_Eqn_Rpt(tauc, tauv, min 1 max 2) 'Parameters Eq.(5-7) local para = Abs((costhetarot*sinphirot)/ax); local parb = Abs((sinthetarot*sinphirot)/bx); local parc = Abs((cosphirot)/cx); 'Parameters chapter 2.4 local pard = Sqrt(para^2+parb^2); local pare = pard/parc; local parf = Sqrt(1-pare^2); 'Here comes the tricky part 'Ellipsoid Eq.(4) local ev = If(modx < 2,CeV(tauc,tauv) * (9/8) * sizeL, 'Cylinder If(modx < 3, 'Eq. (11) If(pard <= 0.0001,CeV(tauc, tauv) * sizeL, 'Eq. (12) If(pard < parc, CeV(tauc, tauv) * ( (2/parc) - (8/(3.14159265358979*pard)) * (-(2/3)+(2/3)*parf+((pare*pare*parf)/3)+(pare*ArcSin(pare)))+ ((2*parc)/(3.14159265358979*pard*pard))*(((pare*parf)/2)-((1/2)*(ArcSin(pare)))+(pare*pare*pare*parf)+(2*pare*pare*ArcSin(pare)))), 'Eq. (13) CeV(tauc, tauv)*((16)/(3*3.14159265358979*pard))-((2*parc)/(4*pard*pard)))), 'Cuboid 'Eq. (8) If(And(para>=parb, para>=parc),CeV(tauc, tauv)*(2/para)*(1-(parb/(3*para))-(parc/(3*para))+((parb*parc)/(6*para*para))), 'Eq. (9) If(And(parb>=para, parb>=parc),CeV(tauc, tauv)*(2/parb)*(1-(para/(3*parb))-(parc/(3*parb))+((para*parc)/(6*parb*parb))), 'Eq. (10) CeV(tauc, tauv)*(2/parc)*(1-(para/(3*parc))-(parb/(3*parc))+((para*parb)/(6*parc*parc)))))) ); 'IB calculations Eq. (1-3) local ibV = 1/ev; local ibL = 1/(2*sizeL); 'ibG is 0 if ibV < ibL local ibG = If (ibV > ibL, ((Sqrt((4*ibV^2)-(4.27679864*ibV*ibL)+(0.27679864*ibL^2)))/2), 0); 'CS_G macro local sizeG = ibG/(Sqrt(3.14159265358979/(4*Ln(2)))); gauss_fwhm = 0.1 Rad Lam sizeG / Cos(Th); 'For Output purposes only local taux = CeV(tauc, tauv); } macro adp_no_limits { move_to u11 min = Val - .05; max = Val + .05; move_to u22 min = Val - .05; max = Val + .05; move_to u33 min = Val - .05; max = Val + .05; move_to u12 min = Val - .05; max = Val + .05; move_to u13 min = Val - .05; max = Val + .05; move_to u23 min = Val - .05; max = Val + .05; } } macro PlotAnisoCS { 'ellipsoid normals_plot = If(modx < 2, sizeL*(3/4), 'cylinder If(modx < 3, 1/( ( ((Abs((1/cx)*cosphirot))^(10))+((Abs((sinphirot))^(10))*( ((Abs((1/bx)*sinthetarot))^(2)) + ((Abs((1/ax)*costhetarot))^(2)) )^(5) ) )^(1/10) ), 'cuboid 1/( ( ((Abs((1/cx)*cosphirot))^(10))+((Abs((sinphirot))^(10))*( ((Abs((1/bx)*sinthetarot))^(10)) + ((Abs((1/ax)*costhetarot))^(10)) ) ) )^(1/10) ) ) ); } macro peak_type_FullsplitPV(consL, consLv, consR, consRv, szL, szLv, szR, szRv, stnL, stnLv, stnR, stnRv, lorL, lorLv, lorR, lorRv, lorL1, lorL1v, lorR1, lorR1v, lorL2, lorL2v, lorR2, lorR2v) { #m_argu consL #m_argu consR #m_argu szL #m_argu szR #m_argu stnL #m_argu stnR #m_argu lorL #m_argu lorR #m_argu lorL1 #m_argu lorR1 #m_argu lorL2 #m_argu lorR2 If_Prm_Eqn_Rpt(consL, consLv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(consR, consRv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(szL, szLv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(szR, szRv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(stnL, stnLv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(stnR, stnRv, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(lorL, lorLv, min 0.0001 max 1, del 0.001) If_Prm_Eqn_Rpt(lorR, lorRv, min 0.0001 max 1, del 0.001) If_Prm_Eqn_Rpt(lorL1, lorL1v, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(lorR1, lorR1v, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(lorL2, lorL2v, min 0.0001 max = 2 Val + .1;, del 0.001) If_Prm_Eqn_Rpt(lorR2, lorR2v, min 0.0001 max = 2 Val + .1;, del 0.001) peak_type spv spv_h1 = CeV(consL,consLv) + CeV(stnL,stnLv) Tan(Th) + CeV(szL,szLv)/Cos(Th) ; spv_h2 = CeV(consR,consRv) + CeV(stnR,stnRv) Tan(Th) + CeV(szR,szRv)/Cos(Th) ; spv_l1 = CeV(lorL,lorLv) + CeV(lorL1,lorL1v) Tan(Th) + CeV(lorL2,lorL2v) / Cos(Th); spv_l2 = CeV(lorR,lorRv) + CeV(lorR1,lorR1v) Tan(Th) + CeV(lorR2,lorR2v) / Cos(Th); } macro Get_Distance(sites, c, d_calc) { #m_ifarg c "" #m_unique_not_refine c #m_endif sites_distance c load site_to_restrain { sites } prm = c; : d_calc } macro Get_Angle(sites, c, ang_calc) { #m_ifarg c "" #m_unique_not_refine c #m_endif sites_geometry c load site_to_restrain { sites } prm = Sites_Geometry_Angle(c); : ang_calc }