;+
;
; Project   : s3drs reconstruction software
;                   
; Name      : (various test routines)
;               
; Purpose   : part of the s3drs test kit
;               
; Explanation: (archived here for completeness, not intended for users use)
;               
; Use       : (self-documented, see notes in code)
;    
; Written     : Sandy Antunes, NRC, 2005-2009
;               
;-            
; Weighting test, works!


nd=2
N=32

; generate some fake data
buildcloud,33,N=N,out=1
timage=readcloudcube(33)
modparam=timage2modparam(timage)
prep_3drv,pxwrapper,views=nd
datum=rtwl_wrapper(pxwrapper,modparam)
tv_multi,datum,title='original data loaded',rescale=128/N

; prepare a pixon object
prep_3drv,pxwrapper,views=nd
p_generate_pixon,'classic',pxwrapper,oPixon
oData = oPixon->get(/oData)
for m=0,nd-1 do oData[m]->Set, data=PTR_NEW(datum[*,*,m]);

; add poisson6 noise
Smax=10
Smin=1
Dfmin=0
Dfmax=max(datum)
a = (Dfmax-Dfmin)/(Smax^2-Smin^2)
Nb = (Dfmax*Smin^2-Dfmin*Smax^2)/(Dfmax-Dfmin)
sigma= sqrt(a*(datum+a*Nb))
oPixon->Replace,sigma=sigma,mxit=100

; add weighting
weight =  make_array(N,N,Nd,value=1);
weight[0:N/2,0:N/2,*]=0; ignore one quarter of each datum
oPixon->Set,wt=weight

; run pseudo
;oldimage=*((oPixon->get(/image))[0])
oldimage=oPixon->get(/image,/nop)
oPixon->pseudo
newimage=oPixon->get(/image,/nop)

; plot results
threeview,oldimage,rescale=128/N,title='input image'
threeview,newimage,rescale=128/N,title='solution image'

; get original data
pdata = oPixon->Get(/data)
data=dblarr(n,n,nd);
for i=0,nd-1 do data[*,*,i]=*pdata[i]
; get the reconstructed images
mdata=dblarr(n,n,nd);
for m=0,nd-1 do mdata[*,*,m] = oData[m]->dsf(newimage)
tv_multi,datum,title='original data refetched',rescale=128/N
tv_multi,mdata,title='solution data',rescale=128/N