hmmm... that dosent post very well. I'd provide a link but i'm
currently without a place to really post the code.
here are the rough examples. at least these may be useful to someone
doing some 4am hacking.
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; ************************************************************ ****
;; code supplied as is. use at own risk, no warranties expressed,
;; implied, or infered.
;; ************************************************************ ***
!p.charsize=1.5
IF 0 THEN BEGIN
x=[0,1,1,0,.2, 1, 0,-.2,1.60001]
y=[0,0,1,1,.2,-1,-1,-.2, 0 ]
v=[0,1,1,1,.1, 3, 3,-.2,1.6 ]
loadct,39,/silent
plot,x,y,color=0,psym=2
loadct,33,/silent
plots,x,y,color=bytscl(v),psym=2,thick=5
loadct,39
radiusbins=.8D
thetabins=2D
r=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabin s,v=v,/
double)
;;now bins are being passed in from above.
r2=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabi ns,v=v,/
double)
phibins=2
r3=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabi ns,phibins=phibins,v=v,/
double)
z=randomu(seed,9)
r4=variogram(x=x,y=y,z=z,radiusbins=radiusbins,thetabins=the tabins,phibins=phibins,v=v,/
double)
ENDIF
IF 1 THEN BEGIN
nside=72L
smoothwidth=15
x=double(floor(findgen(nside,nside)/nside))
y=reform(transpose(x),[nside^2])
x=reform(x,[nside^2])
; z=reform(findgen(nside,nside),[nside^2])
v=randomn(seed,nside,nside)
window,0,xs=1200,ys=800
!p.multi=[0,3,2]
tvscale,v,/nointerp
!p.multi[0]=!p.multi[0]-2
v=reform(v,[nside^2])
radiusbins=5.D ;width
thetabins=3.D ;number of divisions
r=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabin s,cum=cumu,
v=v,/double)
nbins=n_elements(radiusbins)-2
binstarts=radiusbins[0:nbins]
plot,binstarts,r[*,0],yr=[0,max(r)],/ys
FOR d=0,n_elements(thetabins)-2 DO
oplot,binstarts,r[*,d],color=d*254/2., thick=2
!p.multi[0]=!p.multi[0]+3
; plot,cumu[*,0]
; FOR d=0,n_elements(thetabins)-2 DO oplot,cumu[*,d],color=d*50.,
thick=2
;; copy the original data
v0=reform(v,[nside,nside])
v=smooth(v0,floor(smoothwidth*[1.,1./smoothwidth]),/edge) ; smooth
anisotropically
tvscale,reform(v),/nointerp
!p.multi[0]=!p.multi[0]-2
v=reform(v,[nside^2])
r=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabin s,cum=cumu,
v=v,/double)
nbins=n_elements(radiusbins)-2
binstarts=radiusbins[0:nbins]
plot,binstarts,r[*,0],yr=[0,max(r)],/ys
FOR d=0,n_elements(thetabins)-2 DO
oplot,binstarts,r[*,d],color=d*254/2., thick=2
!p.multi[0]=!p.multi[0]+3
v=smooth(v0,floor(smoothwidth),/edge) ;isotropic
tvscale,reform(v),/nointerp
!p.multi[0]=!p.multi[0]-2
v=reform(v,[nside^2])
r=variogram(x=x,y=y,radiusbins=radiusbins,thetabins=thetabin s,cum=cumu,
v=v,/double)
nbins=n_elements(radiusbins)-2
binstarts=radiusbins[0:nbins]
plot,binstarts,r[*,0],yr=[0,max(r)],/ys
FOR d=0,n_elements(thetabins)-2 DO
oplot,binstarts,r[*,d],color=d*254/2., thick=2
ENDIF
IF 0 THEN BEGIN
nside=50^2
x=randomu(seed,nside)-.5
y=randomu(seed,nside)-.5
z=randomu(seed,nside)-.5
v=sqrt( x^2. + y^2. + z^2. )
; x=x[sort(v)]
; y=y[sort(v)]
;; z=z[sort(v)]
; v=v[sort(v)]
loadct,39,/silent
window,1,xs=800,ys=800
SURFACE, DIST(5), /NODATA, /SAVE, XRANGE=[0,1]-.5, $ ;thanks David,
the RSI examples are worthless.
YRANGE=[0,1]-.5, ZRANGE=[0, 1]-.5, XSTYLE=1, $
YSTYLE=1, ZSTYLE=1, CHARSIZE=1.5, $
POSITION=[0.1, 0.1, 0.95, 0.95, 0.1, 0.95], $
XTICKLEN=1, YTICKLEN=1, XGRIDSTYLE=1, YGRIDSTYLE=1
AXIS, XAXIS=1, /T3D, CHARSIZE=1.5
AXIS, YAXIS=1, /T3D, CHARSIZE=1.5
phi = Findgen(32) * (!PI * 2 / 32.)
phi = [ phi, phi(0) ]
UserSym, Cos(phi), Sin(phi), /Fill
loadct,33,/silent
PLOTS, x, y, z, PSYM=8, COLOR=bytscl(v), SYMSIZE=2.5, /T3D
plots, [0,0],[0,0],[-1,1]*.5,/t3d,thick=2
plots, [0,0],[-1,1]*.5,[0,0],/t3d,thick=2
plots, [-1,1]*.5,[0,0],[0,0],/t3d,thick=2
; FOR pp=0,nside-1 DO plots,[0,x[pp]],[0,y[pp]],[0,z[pp]],/t3d
loadct,39,/silent
radiusbins=.05
thetabins=[-1D*!dpi/2.D, .1D, !dpi/2.D] ;2 ;trying regular and
irregularly specified bins
phibins=[0.D, (!dpi/(2.D))-1D, !dpi] ;2
r=variogram(x=x,y=y,z=z,$
radiusbins=radiusbins,thetabins=thetabins,phibins=phibins,$
cum=cumu, v=v,/double)
nbins=n_elements(radiusbins)-2
binstarts=radiusbins[0:nbins]
window,2,xs=800,ys=800
plot,binstarts,r[*,0,0],yr=[0,max(r)]
FOR t=0,1 DO FOR p=0,1 DO
oplot,binstarts,r[*,t,p],color=((t*2)+p)*254/3.
ENDIF
END
|