-
Notifications
You must be signed in to change notification settings - Fork 2
/
restore.f90
80 lines (62 loc) · 2.05 KB
/
restore.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
subroutine restore_decomp (img,outimg,nx,ny,oimg,work,wksp,&
scalemin,scalemax,savehist,shdescriptor,qnoise,nimg,outkey)
implicit none
integer :: nx,ny
real, dimension(nx*ny) :: img, outimg, wksp
integer :: scalemin,scalemax
integer (kind=1), dimension(nx*ny,*) :: savehist
logical :: qnoise
real, dimension(nx*ny,*) :: nimg
character(len=20), dimension(300) :: shdescriptor
integer :: nrec
real, dimension(nx*ny) :: oimg,work
integer :: w1,w2
integer :: irec,iscale,iiter,i, iscaleold
integer :: kernel
character (len=*) :: outkey
character (len=20) :: affix
character (len=200) :: outname
!c ///
read(shdescriptor(1),*)w1,w2,nrec, kernel
call msarith (oimg,'=',0.0,nx*ny)
!c ///
iscaleold = -1
do irec=1,nrec
read (shdescriptor(irec+1),*) iscale,iiter
if (iscale.ne.iscaleold.and. &
(iscale.ge.scalemin.and.iscale.le.scalemax)) then
if (iscaleold.gt.0) then
affix=' '
write (affix,'(i2)') iscaleold
call rmblanks(affix)
outname = outkey
call strcat (outname,'.')
call strcat (outname,affix)
call write_fits_image (outname,outimg,nx,ny,'e','e')
endif
iscaleold=iscale
call msarith (outimg,'=',0.0,nx*ny)
endif
print*,'Restoring: ',iscale,iiter
call marith (work,'=',img,'-',oimg,nx*ny)
call conv_wavelet (kernel,work,work,wksp,nx,ny,iscale,.true.)
call marith (work,"=",work,"-",wksp,nx*ny)
do i=1,nx*ny
if (savehist(i,irec).eq.0) work(i)=0.0
enddo
call marith (oimg,'=',oimg,'+',work,nx*ny)
if (iscale.ge.scalemin.and.iscale.le.scalemax) then
call marith (outimg,'=',outimg,'+',work,nx*ny)
endif
enddo
if (iscaleold.gt.0) then
affix=' '
write (affix,'(i2)') iscaleold
call rmblanks(affix)
outname = outkey
call strcat (outname,'.')
call strcat (outname,affix)
call write_fits_image (outname,outimg,nx,ny,'e','e')
endif
return
end subroutine restore_decomp