; .run getgrid first
;
; plot histograms of the grid
; wedge data should be in dar
if 0  then begin
histMinCnt=5
stpAr=[.3,.6,.9,1.2]
nstp=n_elements(stpar)
nscans=9
sclp=1000. ; to mm
minR=0.
maxR=4.
minA=-4.
maxA=4.
bin=.2
endif
;goto,plotit1
;goto,avgscans
;goto,plotit
start=1
;
; for scan 10 get the raw errors.
;
plotit:
;
	if hard then pscol,'wedge_gridhist.ps',/full
;
; 	show histogram for scan 10 grid
;
	;
; for scan 10 get the raw errors.
;
	cs=1.4
	font=1
	csn=1.5
	wedge_getgrid,wi,dar[0],.3,gdI,xyz=xyz
	Rm=reform(blmcmpradius(xyz,center=wI.center))
	Rerr= Rm - wi.r0
	bin=.1
	max=10.
	min=-10
	h=histogram(rerr*1000.,bin=bin,max=max,min=min,loc=hloc)
	nterms=4
	a=dblarr(nterms)
	a[0]=max(h)
	a[1]=mean(h)
	aa=rms(rerr*1000,/quiet)
	a[2]=aa[1]
	a[3]=0.
	yfit=gaussfit(hloc,h,a,nterms=nterms,yerror=yerror)
	xtit='Radial error [mm]'
	ytit='counts'
	tit='11mar20 histogram of radial err over scan 10 grid'
	!p.multi=0
	ver,0,1100
	hor,min,max
	plot,hloc,h,$
		chars=cs,font=font, $
  		xtit=xtit,ytit=ytit,$
        tit=tit
	oplot,hloc,yfit,col=colph[2]
	ln=3
	xp=.04
	csn=1.5
	scl=.8
note         ,ln,'gaussfit to histogram',chars=csn,font=font,col=colph[2],xp=xp
	lin=string(format='("mean = ",f5.2," mm")',a[1])
note,ln+1*scl,lin,chars=csn,font=font,col=colph[2],xp=xp
	lin=string(format='("sigma= ",f5.2," mm")',a[2])
note,ln+2*scl,lin,chars=csn,font=font,col=colph[2],xp=xp

; 	plot mean rms for entire grid  by scan
	!p.multi=[0,1,2]
	cs=1.4
	x=findgen(nscans)+10 
	hor
	ver
	sym=-1
	ver,-3,3
	hor,9,19
	xtit='scan'
	ytit='radial err [mm]'
	tit='11mar20 radial err over grid (by scan)'
	plot,x,har[0,*].avgall*sclp,psym=sym,$
			chars=cs,font=font,$
			xtit=xtit,ytit=ytit,$
			tit=tit
    ver,0,4
	tit='rms of Radial err over grid (by scan)'
	plot,x,har[0,*].rmsall*sclp,psym=sym,$
			chars=cs,font=font,$
			xtit=xtit,ytit=ytit,$
			tit=tit
;
; 	now look at rms by grid size for each 

ver,1,800
hor,0,4	 	
	!p.multi=[0,1,3]
	cs=2.2
	ytit='counts'
	xtit='rms [mm]'
	tit='histogram of Radial err rms computed by grid bin size'
	for iscan=0,nscans-1  do begin &$
		lscan=string(format='(" Scan:",i2)',har[0,iscan].scan ) &$
		plot,[0,1],[0,1],/nodata,/ylog,$
			chars=cs,font=font,$
			xtit=xtit,ytit=ytit,$
			tit=tit + lscan &$
		for istp=0,nstp-1 do begin &$
		   icol=istp+1 &$
			oplot,hlocr,har[istp,iscan].hrms,psym=sym,col=colph[icol] &$
		endfor &$
		if iscan mod 3 eq 0 then begin
			ln=2.5
			xp=.6
		    nscl=.7	
			for j=0,nstp-1 do begin
				note,ln + j*nscl,string(format='("binSize:",f3.1," m")',stpar[j]),$
					chars=csn,font=font,xp=xp,col=colph[j+1]
			endfor
		endif
	endfor
	if hard then hardcopy
	x
	ldcolph
;
; for a given scan:
; 1. start with the grid of averaged points
; 2. remove the mean value of scan from entire grid.
; 3. compute the rms between scans 
; 4. do this with and without removing the means;
;
avgScans:
stp=.3
	for iscan=0,nscans-1 do begin &$
    	print,"start scan:",iscan &$
        wedge_getgrid,wi,dar[iscan],stp,gdI &$
		if iscan eq 0 then begin 
			a={ npnts: gdI.npnts,$
				nx   : gdI.nx,$
	            ny   : gdI.ny,$
               avgAll: gdI.avgall,$
               rmsAll: gdI.rmsall,$
			   xar   : gdI.xar,$
			   yar   : gdI.yar,$
              cntAr  : gdI.cntAr,$
              rmsAr  : gdI.rmsAr,$
			  avgAr  : gdI.avgAr}
			gdIAr=replicate(a,nscans)
			gdI1=gdiAr[0]
		endif
		struct_assign,gdI,gdI1
		gdiAr[iscan]=gdi1
		print,rms(gdiar[iscan].avgAr)
 	endfor
	stop
;
;	now remove the mean value from each bin err
;
	for i=0,nscans-1 do print,rms(gdiar[i].avgAr - mean(gdiar[i].avgar)) 

	avgAr1Avg=dblarr(nx,ny)
	avgAr1Med=dblarr(nx,ny)
	rmsAr1Avg=dblarr(nx,ny)
	rmsAr1Med=dblarr(nx,ny)
	avgAr=dblarr(nx,ny,nscans)
	medAr=dblarr(nx,ny,nscans)
	for i=0,nscans -1 do begin &$
		avgAr[*,*,i]=gdiar[i].avgar - mean(gdiar[i].avgAr) &$
		medAr[*,*,i]=gdiar[i].avgar - median(gdiar[i].avgAr) &$
	endfor
	for ix=0,gdi1.nx-1 do  begin &$
		for iy=0,gdi1.ny-1  do begin &$
			a=rms(avgAr[ix,iy,*],/quiet) &$
			avgAr1Avg[ix,iy]=a[0] &$
			rmsAr1Avg[ix,iy]=a[1] &$
			a=rms(medAr[ix,iy,*],/quiet) &$
			avgAr1Med[ix,iy]=a[0] &$
			rmsAr1Med[ix,iy]=a[1] &$
		endfor &$
	endfor
;
		harNdc=lonarr(40,nscans)
		for i=0,nscans-1 do begin &$
        	hArNDc[*,i]=histogram((avgar[*,*,i] - mean(gdiar[i].avgAr))*sclp,$
			bin=bin,min=minA,max=maxA,loc=hlocE) &$
		endfor
		
	    binA=.05
        havg=histogram(avgar1avg*sclp,bin=bin,min=minA,max=maxA,loc=hlocE)
        hrms=histogram(rmsar1avg*sclp,bin=binA,min=minR,max=maxR,loc=hlocRA)
plotit1:
	if hard then pscol,'wedge_gridavg.ps',/full
!p.multi= 0
	cs=1.2
	font=1
	csn=1.5
	psym=-1
	hor,minr,maxr
	ver,0,800
	xtit='rms [mm]'
	ytit='counts'
	tit='11mar20 grids. histogram of  rms errors within .3m bins'
	stripsxy,hlocr,har[0,*].hrms,0,0,psym=sym,$
		xtit=xtit,ytit=ytit,tit=tit,chars=cs,font=font
	oplot,hlocRA,hrms,psym=sym,col=colph[2]
	fl1=1.199392
	fl2=0.39490
	flag,fl1,col=colph[3]
	flag,fl2,col=colph[3]
	eps=.02
	xyouts,fl2+eps,500,'.395mm',chars=csn,font=font,col=colph[3]
	xyouts,fl1+eps,100,'1.200mm',chars=csn,font=font,col=colph[3]
	
	ln=2
	xp=.3
	scl=.7
	note,ln,'9 scans (no averaging) (about 30pnts/bin)',chars=csn,font=font,$
			xp=xp
	note,ln+1*scl,"avg each bin, Then compute rms between scans",$
		chars=csn,font=font,xp=xp,col=colph[2]
	lin="rms improvement: 1.2mm/.395mm=3.0"
	ln=5
	xp=.4
	note,ln,lin,chars=csn,font=font,xp=xp
	note,ln+1,"expected rms improvement: sqrt(30)=5.5",chars=csn,font=font,xp=xp

	if hard then hardcopy
	x
	ldcolph
end