4Misc_StartP#4Platform@9VersionCheck xHH@Rg(HHdh xHH@Rg(HHdh x HH@Rg(HHdh ^Graph*xpWDashSettings#  !4 4 4 4 4 4 whomeCCdMacintosh HD:Users:lingyuan:Desktop:mRecentWindows9Auxiliaryd3dErrors.ihfHaldane ModelQi-Wu-Zhang Model 4Misc_EndP#4XOPState_StartP#4XOPState_EndP#\,A?E'm1:!J*_d TXET0RGI(_Qi-Wu-Zhang ModelQi-Wu-Zhang Model#pragma TextEncoding = "UTF-8" #pragma rtGlobals=3 // Use modern global access method and strict wave access #pragma DefaultTab={3,20,4} // Set default tab width in Igor Pro 9 and later /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// //**************************** Qi-Wu-Zhang Model ****************************************** /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////////////////// //X.-L. Qi, Y.-S. Wu, and S.-C. Zhang, Topological quantization of the spin Hall effect in two-dimensional paramagnetic semiconductors, Phys. Rev. B 74, 085308 (2006). //Ref: https://journals.aps.org/prb/abstract/10.1103/PhysRevB.74.085308 Function/Wave QWZ(A,B,m,kx,ky,mu) variable A,B,m,kx,ky,mu //** Pointer String for called Functions string ps = "ps_QWZ" string result_C = "result_C" //** Create QHZ Hamiltonian variable pf1,pf2,pf3,pf4,pf5,pf6 pf1 = A*sin(kx) //From A*kx; kx --> sin(kx) pf2 = A*sin(ky) pf3 = m+B*(2-cos(kx)-cos(ky)) //Tranform from m+B*(kx^2+ky^2); kx^2 --> 1-cos(kx) pf4 = -mu wC("pf1_QWZ",{{cmplx(pf1,0)}}) wC("pf2_QWZ",{{cmplx(pf2,0)}}) wC("pf3_QWZ",{{cmplx(pf3,0)}}) wC("pf4_QWZ",{{cmplx(pf4,0)}}) Wn("ps_QWZ",{{1},{2},{3},{0}}) automatrixC("pf1_QWZ;pf2_QWZ;pf3_QWZ;pf4_QWZ",$ps) wave/C result_Cw = $result_C duplicate/C/o result_Cw QWZH killwaves $result_C $ps return QWZH //** Error message, no stop but Continuing execution String msg Variable err msg=GetErrMessage(GetRTError(0),3); err=GetRTError(1) if (err != 0) //Print "Error in Demo: " + msg //Print "Continuing execution" endif end Function QWZSpectralFunction(A,B,m,mu,selfimg,numxy) variable A,B,m,mu variable selfimg variable numxy //= 100 variable/G zn_cons = 2 //** Make void matrix for spectral function make/o/n=(numxy,numxy,500) modelresult_QWZ setscale/i x,-pi,pi,"",modelresult_QWZ setscale/i y,-pi,pi,"",modelresult_QWZ setscale/i z,-3,3,"",modelresult_QWZ modelresult_QWZ = 0 variable roff = dimoffset(modelresult_QWZ,2) variable rdelta = dimdelta(modelresult_QWZ,2) //** Variable for Spectral Function variable i,j,selfreal,t selfreal = 0 variable kx variable ky //** Pointer String for called Functions string W_eigenvalues="W_eigenvalues" string QWZH = "QWZH" //** Start loop for (kx, ky) i=0 do j=0 do //** Define kx, ky by void matrix shape kx = dimoffset(modelresult_QWZ,0)+i*dimdelta(modelresult_QWZ,0) ky = dimoffset(modelresult_QWZ,1)+j*dimdelta(modelresult_QWZ,1) QWZ(A,B,m,kx,ky,mu) wave/C QWZHw = $QWZH //** Solve eigenvalue of QWZ Hamiltonian matrixEigenV QWZHw wave/C n=$W_eigenvalues make/n=(dimsize(QWZHw,0))/o sorteigen sorteigen[0]= real(N[0]) sorteigen[1]= real(N[1]) //Make 3D Spectral Function t=0 do modelresult_QWZ[i][j][]+=(-1/pi)*selfimg/(selfimg^2+((roff+r*rdelta)-selfreal-sorteigen[t])^2) t+=1 while (t < dimsize(sorteigen,0)) j+=1 while(j < dimsize(modelresult_QWZ,1)) i+=1 while(i < dimsize(modelresult_QWZ,0) ) //** Call smart 3D displayer duplicate/o modelresult_QWZ Z2_001 //killwaves modelresult_QWZ sorteigen QWZHw n execute "d3dsimu(\"Z2_001\",2)" //** Modify the 3D player for this special simulation use (band structure) variable/G divcolor_cons divcolor_cons = 1 wavestats/Q $tpw() variable rangls = max(abs(V_max),abs(V_min)) ModifyImage/W=$grabwinchild(tpw()) $tpw() ctab= {-rangls,rangls,RedWhiteBlue,0}; ColorScale/W=$grabwinchild(tpw())/C/N=textcc/X=-21.00/Y=5.00/F=0 heightPct=30,nticks=5,minor=1,tickLen=1.00,image=$tpw()//;ColorScale/W=$grabwinchild(Fexyd)/C/N=text0/F=0/A=MC/Y=-120 image=$Fexyd variable/G normornot_3dplot =2 Button launchLinecut title="Linecut",size={65,14},fSize=10,proc=ButtonProc_Cons3dplotZ end Function QWZBerryCurvature(A,B,m,mu,numxy) variable A,B,m,mu variable numxy //= 100 variable kx variable ky variable i,j Variable timerRefNum timerRefNum = StartMSTimer //** Make void matrix for spectral function make/o/n=(numxy,numxy,1) voidberry_QWZ setscale/i x,-pi,pi,"",voidberry_QWZ setscale/i y,-pi,pi,"",voidberry_QWZ //setscale/i z,-2.5,4.5,"",voidberry_QWZ voidberry_QWZ = 0 //** Make void matrix for Berry curvature band #1 (lower) and #2 (higher) make/o/n=(numxy,numxy) BC_QWZ1 setscale/i x,-pi,pi,"",BC_QWZ1 setscale/i y,-pi,pi,"",BC_QWZ1 BC_QWZ1 = 0//cmplx(0,0) make/o/n=(numxy,numxy) BC_QWZ2 setscale/i x,-pi,pi,"",BC_QWZ2 setscale/i y,-pi,pi,"",BC_QWZ2 BC_QWZ2 = 0//cmplx(0,0) //** Pointer String for called Functions string W_eigenvalues="W_eigenvalues" string QWZH = "QWZH" string M_R_eigenVectors = "M_R_eigenVectors" string M_product = "M_product" //** variable for Berry curvature variable E1k, E2K //1 for low band, 2 for high band variable delta = 10^-6 variable/C Berryterm1A = cmplx(0,0) variable/C Berryterm1B = cmplx(0,0) variable/C Berryterm2A = cmplx(0,0) variable/C Berryterm2B = cmplx(0,0) //** Start loop for (kx, ky) i=0 do j=0 do //** Define kx, ky by void matrix shape kx = dimoffset(voidberry_QWZ,0)+i*dimdelta(voidberry_QWZ,0) ky = dimoffset(voidberry_QWZ,1)+j*dimdelta(voidberry_QWZ,1) QWZ(A,B,m,kx,ky,mu) wave/C QWZHw = $QWZH //** Solve eigenvalue of QWZ Hamiltonian matrixEigenV/R QWZHw wave/C n=$W_eigenvalues wave/C eignevector=$M_R_eigenVectors make/n=(dimsize(QWZHw,0))/o sorteigen sorteigen[0]= real(N[0]) sorteigen[1]= real(N[1]) //**Sort eigenvalue and eigenvector make/o/n=(dimsize(QWZHw,0)) sortindex sortindex=p sort sorteigen sortindex //the value indicates the column number in M_R_eigenVectors belong to each eigenvalue labeled in sorteigen sort sorteigen sorteigen //** Calculat Berry Curvature //## E1,k = sorteigen[0] //## |u1,k> = eignevector[][sortindex[0]] //## E2,k = sorteigen[1] //## |u2,k> = eignevector[][sortindex[1]] //** get eigenvalue E1k = sorteigen[0] E2k = sorteigen[1] //** make |u1,k> //|ket> make/n=(dimsize(QWZHw,0),1)/o/C u1k_ket u1k_ket[] = eignevector[p][sortindex[0]] // //|ket> make/n=(dimsize(QWZHw,0),1)/o/C u2k_ket u2k_ket[] = eignevector[p][sortindex[1]] // = eignevector[][sortindex[0]] //## E2,k = sorteigen[1] //## |u2,k> = eignevector[][sortindex[1]] //** get eigenvalue E1k = sorteigen[0] E2k = sorteigen[1] //** make |u1,k> //|ket> make/n=(dimsize(QWZHw,0),1)/o/C u1k_ket u1k_ket[] = eignevector[p][sortindex[0]] // //|ket> make/n=(dimsize(QWZHw,0),1)/o/C u2k_ket u2k_ket[] = eignevector[p][sortindex[1]] //V_min) lc = W_coefw[2]-0.5*tt*sigma else lc =V_min endif if (W_coefw[2]+0.5*tt*sigma < V_max) lh = W_coefw[2]+0.5*tt*sigma else lh =V_max endif ModifyImage/W=$grabwinnonew(nameofwave(name)) $nameofwave(name) ctab= {lc,lh,VioletOrangeYellow,0} //PRINT num2str(W_coefw[2]-1.5*sigma)+" "+ num2str(W_coefw[2]+1.5*sigma) //print num2str(W_coefw[2]) //display namehistw end Function gethistgram_npcolor(name) string name string histn = "Hist_"+name make/o/N=1000 $histn wavestats/Q $name variable Vmax = V_max Variable Vmin = V_min variable mdian = median($name) Histogram/B={(mdian-(Vmax-Vmin)/2),(V_max-V_min)/1000,1000} $name,$histn duplicate/o $histn fithttemp //CurveFit/Q gauss, $histn/D=fithttemp variable V_FitError = 0 CurveFit/Q gauss, $histn/D=fithttemp if (V_FitError != 0) endif killwaves fithttemp string W_coef="W_coef" wave W_coefw = $W_coef variable coverratio //(0,1] coverratio = (2*5*W_coefw[3])/(Vmax-Vmin) Histogram/B={(mdian-coverratio*(Vmax-Vmin)/2),coverratio*(V_max-V_min)/200,200} $name,$histn //CurveFit/Q gauss, $histn/D//=fithttemp V_FitError = 0 CurveFit/Q gauss, $histn/D//=fithttemp if (V_FitError != 0) endif end Function/s grabwinnonew(name) string name string fulllist = WinList("*", ";","WIN:1") string nn,waveong,cmdn,out out = "" variable i for(i=0; iV_min) lc = W_coefw[2]-0.5*tt*sigma else lc =V_min endif if (W_coefw[2]+0.5*tt*sigma < V_max) lh = W_coefw[2]+0.5*tt*sigma else lh =V_max endif ModifyImage/W=$grabwinnonew(nameofwave(name)) $nameofwave(name) ctab= {lc,lh,VioletOrangeYellow,1} //PRINT num2str(W_coefw[2]-1.5*sigma)+" "+ num2str(W_coefw[2]+1.5*sigma) //print num2str(W_coefw[2]) //display namehistw end Function/s grabwin(name) string name string fulllist = WinList("*", ";","WIN:1") string nn,waveong,cmdn,out out = "" variable i for(i=0; i= 0 && round(qq) <= Ny-1 ) //** condition of Q in range pps = i break endif i+=1 while(i < Nx) return pps end //***************************************************** //## 03 //** Determine the Last P for assignment loop Function findendppf(mat,angle,Zn,addY) String mat variable angle variable Zn variable addY //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //print num2str(xn)+" "+num2str(yn) //** Initial layer dimsize wave matw = $mat variable Nx,Ny,Nz Nx = dimsize(matw,xn) Ny = dimsize(matw,yn) Nz = dimsize(matw,zn) //** Search for the End P variable pp, qq,ppe variable i ppe=nan i=Nx-1 //## Search from the largest P, calculate when Q is in range, that position is ppe do qq = tan(angle*pi/180)*(i-round(Nx/2))+round(Ny/2)+addY //[Formula of Rotated HLine] //print qq if(round(qq) <= Ny-1 && round(qq) >= 0) //** condition of Q in range ppe = i break endif i-=1 while(i >= 0) return ppe end //***************************************************** //## 04 //** Find the addY Range for rotated H-Linecut Function/Wave findrangeforangle_LHf(mat,angle,Zn) String mat variable angle variable Zn //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //** Search for high limit variable addY, addYmax, addYmin,body addY =0 do body = mod(Round(findendppf(mat,angle,Zn,addY)),1) if (body != 0) addYmax = addY-1 break endif addY+=1 while(addY < 3*dimsize($mat,yn)) //** Search for low limit addY =0 do body = mod(Round(findendppf(mat,angle,Zn,addY)),1) if (body != 0) addYmin = addY+1 break endif addY-=1 while(addY > -3*dimsize($mat,yn)) //** Return the limit wave make/N=2/o Hcut_Vrange //print addYmin //print addYmax Hcut_Vrange={addYmin,addYmax} return Hcut_Vrange end //////////////////////////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////////////////////////// //** INTERACTING FUNCTIONAL: Extract rotated Vertical Linecut //////////////////////////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////////////////////////// //%%%%%%%%%%%%** get rotated V-linecut, addX here is counted from center. //%%%%%%%%%%%%** V-linecut can move in X direction, addX is DP(Int) //////////////////////////////////////////////////////////////////////////////////////////////////////////// //***************************************************** //## 01 //** Main Function Function anglelinecutVf(mat,angle,Zn,addX,normornot,smornot) String mat variable angle variable Zn variable addX variable normornot variable smornot //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //print num2str(xn)+" "+num2str(yn) //** Initial layer dimsize wave matw = $mat variable Nx,Ny,Nz Nx = dimsize(matw,xn) Ny = dimsize(matw,yn) Nz = dimsize(matw,zn) variable pp, qq //** Determine Limit for Assignment Loop //## Play with vertical Line //## Logistically, we swap P and Q //## that mean, set Q as independent variable //## and check P the dependent variable //** Determine the First and Last Q for value assignment loop variable qqs,qqe,pps,ppe,xxs,xxe,yys,yye //pp = tan(angle*pi/180)*(qq-round(Ny/2))+round(Nx/2)+addX qqs=findstartqqf(mat,angle,Zn,addX) //First Q qqe=findendqqf(mat,angle,Zn,addX) //Last Q //** Calculate P and X, Y yys = dimoffset(matw,yn)+qqs*dimdelta(matw,yn) yye = dimoffset(matw,yn)+qqe*dimdelta(matw,yn) pps = tan(-angle*pi/180)*(qqs-round(Ny/2))+round(Nx/2)+addX ppe = tan(-angle*pi/180)*(qqe-round(Ny/2))+round(Nx/2)+addX xxs = dimoffset(matw,xn)+pps*dimdelta(matw,xn) xxe = dimoffset(matw,xn)+ppe*dimdelta(matw,xn) //** Calculate linecut scale variable len = sqrt((xxe-xxs)^2+(yye-yys)^2) string linecutV = "LV_"+mat make/N=((qqe-qqs+1),Nz)/o $linecutV wave linecutVw = $linecutV setscale/p y,dimoffset(matw,zn),dimdelta(matw,zn),"",linecutVw setscale/i x,yys,yys+len,"",linecutVw //** Extract value from 3D matrix //** The value assignment loop is runing along Q variable i i=0 qq=qqs do //qq= round(tan(-angle*pi/180)*(pp-round(Nx/2))+round(Ny/2)+addY) pp = round(tan(-angle*pi/180)*(qq-round(Ny/2))+round(Nx/2)+addX) //[Formula of Rotated VLine] //linecutHw[][i] = matw[pp][qq][p] if(zn == 0) linecutVw[i][] = matw[q][pp][qq] endif if(zn == 1) linecutVw[i][] = matw[pp][q][qq] endif if(zn == 2) linecutVw[i][] = matw[pp][qq][q] endif qq+=1 i+=1 while(qq < qqe+1) //di(linecutHw) //** Shift the zero for FFT linecut string slicename = "Zslice_"+mat wave slicenamew = $slicename variable kk = (xxs-xxe)/(yys-yye) variable bb = xxs-kk*yys variable zx,zy zy = -(bb*kk)/(kk^2+1) zx = bb/(kk^2+1) variable zp,zq zp = round((zx-dimoffset(slicenamew,0))/dimdelta(slicenamew,0)) zq = round((zy-dimoffset(slicenamew,1))/dimdelta(slicenamew,1)) variable dimoff = dimoffset(linecutVw,0) variable dimdel = dimdelta(linecutVw,0) variable correcteddimoffset = -(zq-qqs)*dimdel setscale/p x,correcteddimoffset,dimdel,"",linecutVw //** Shift the zero for FFT linecut //string slicename = "Zslice_"+mat //wave slicenamew = $slicename //variable kk = (xxs-xxe)/(yys-yye) //variable bb = xxs-kk*yys //variable zx,zy //zy = -(bb*kk)/(kk^2+1) //zx = bb/(kk^2+1) //variable zp,zq //zp = round((zx-dimoffset(slicenamew,0))/dimdelta(slicenamew,0)) //zq = round((zy-dimoffset(slicenamew,1))/dimdelta(slicenamew,1)) //correcteddimoffset + dimdelta(linecutHw,0)*zp = 0 //variable correcteddimoffset = -dimdelta(slicenamew,1)*zq //setscale/p x,correcteddimoffset,dimdelta(slicenamew,1),"",linecutVw end //***************************************************** //## 02 //** Determine the First Q for assignment loop Function findstartqqf(mat,angle,Zn,addX) String mat variable angle variable Zn variable addX //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //print num2str(xn)+" "+num2str(yn) //** Initial layer dimsize wave matw = $mat variable Nx,Ny,Nz Nx = dimsize(matw,xn) Ny = dimsize(matw,yn) Nz = dimsize(matw,zn) //** Search for the Start Q variable pp, qq,qqs variable i qqs=nan i=0 //## Search from the smallest Q [i.e. 0], calculate when P is in range, that position is pps do pp = tan(-angle*pi/180)*(i-round(Ny/2))+round(Nx/2)+addX //[Formula of Rotated VLine] if(round(pp) >= 0 && round(pp) <= Nx-1 ) //Condition for P in range qqs = i break endif i+=1 while(i < Ny) return qqs end //***************************************************** //## 03 //** Determine the Last Q for assignment loop Function findendqqf(mat,angle,Zn,addX) String mat variable angle variable Zn variable addX //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //print num2str(xn)+" "+num2str(yn) //** Initial layer dimsize wave matw = $mat variable Nx,Ny,Nz Nx = dimsize(matw,xn) Ny = dimsize(matw,yn) Nz = dimsize(matw,zn) //** Search for the End Q variable pp, qq,qqe variable i qqe=nan i=Nx-1 //## Search from the largest Q, calculate when P is in range, that position is pps do pp = tan(-angle*pi/180)*(i-round(Ny/2))+round(Nx/2)+addX //[Formula of Rotated VLine] //print qq if(round(pp) >= 0 && round(pp) <= Nx-1) //Condition for P in range qqe = i break endif i-=1 while(i >= 0) return qqe end //***************************************************** //## 04 //** Find the addX Range for rotated V-Linecut Function/Wave findrangeforangle_LVf(mat,angle,Zn) String mat variable angle variable Zn //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum //** High Limit variable addX, addXmax, addXmin,body addX =0 do body = mod(Round(findendqqf(mat,angle,Zn,addX)),1) if (body != 0) addXmax = addX-1 break endif addX+=1 while(addX < 3*dimsize($mat,xn)) //** Low Limit addX =0 do body = mod(Round(findendppf(mat,angle,Zn,addX)),1) if (body != 0) addXmin = addX+1 break endif addX-=1 while(addX > -3*dimsize($mat,xn)) //** Return Limit make/N=2/o Vcut_Vrange //print addYmin //print addYmax Vcut_Vrange={addXmin,addXmax} return Vcut_Vrange end //////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////// //** Auxiliary Functions //////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////// Function sumoned(dest,orin,zn) string dest //Output name string orin //Input 3D matrix variable zn //Index of your interesting dimension that keep in the end. //** Auto ordering layer index make/o/N=3 ordernum ={0,1,2} ordernum=abs(ordernum-zn) wavestats/Q ordernum ordernum={0,1,2} DeletePoints V_minloc,1, ordernum wavestats/Q ordernum variable xn = ordernum[V_minloc] variable yn = ordernum[V_maxloc] killwaves ordernum wave orinw = $orin make/N=(dimsize(orinw,zn))/o $dest wave destw = $dest destw=0 variable i,j,k,addk,countnum=0 variable is,ie,js,je variable/G FFTmarque_3dplot// 1 for whole, 2 for in marque if (FFTmarque_3dplot == 1) is = 0 ie = dimsize(orinw,xn) js = 0 je = dimsize(orinw,yn) endif if (FFTmarque_3dplot == 2) getmarquee/W=$winname(0,1) left, bottom variable p1,p2,q1,q2 is = round((v_left-dimoffset($tpw(),0))/dimdelta($tpw(),0)) ie = round((v_right-dimoffset($tpw(),0))/dimdelta($tpw(),0))+1 js = round((v_bottom-dimoffset($tpw(),1))/dimdelta($tpw(),1)) je = round((v_top-dimoffset($tpw(),1))/dimdelta($tpw(),1))+1 //print is, ie ,js, je //print tpw() endif //if (mod(Round(mean(orinw)),1)!=0) i=is//round(dimsize(orinw,xn)/6) do j=js//round(dimsize(orinw,yn)/6) do if(zn== 0) destw[]+=orinw[p][i][j] endif if(zn== 1) destw[]+=orinw[i][p][j] endif if(zn== 2) //destw[]+=orinw[i][j][p] k=0 addk=0 do addk+=orinw[i][j][k] k+=1 while(k jumpjudge) //print i InsertPoints dimsize(jp_p1d,0),1, jp_p1d //print dimsize(jp_p1d,0)-1 jp_p1d[dimsize(jp_p1d,0)-1] = (i+1)*sign(deltawavew[i]) endif i+=1 while (i