generated from erosita/uds
220 lines
6.9 KiB
Fortran
220 lines
6.9 KiB
Fortran
!c===============================================================
|
|
subroutine read_cobe(dglon,dglat,dflux_49,tnrows)
|
|
!c--------------------------------------------------------------
|
|
!c write image to FITS file
|
|
!c
|
|
!c Created: Mon May 16 14:21:02 EDT 1994
|
|
!c--------------------------------------------------------------
|
|
implicit none
|
|
character filename*40
|
|
integer nxm,nx,ny,erange,n,m
|
|
logical anyf
|
|
logical exact,anynull
|
|
|
|
integer lnblnk
|
|
integer bitpix,pcount
|
|
real*8 norma0,norma1,slope0,slope1,dist,e1,e2,texp
|
|
real*8 crval1,crval2,crpix1,crpix2,cdelt1,cdelt2,crota2
|
|
integer naxes(3),naxis,gcount,extend
|
|
integer unit,status,block,extver,hdutype,simple
|
|
character*400 comment
|
|
integer nmapm,nmapx,nmapy
|
|
parameter(nmapm=2000)
|
|
integer nhdu,tnrows,tncols
|
|
|
|
real skyexp(nmapm,nmapm)
|
|
real skydld(nmapm,nmapm)
|
|
real skyerr(nmapm,nmapm)
|
|
real image(nmapm,nmapm)
|
|
|
|
! cobe
|
|
integer TROWMAX,ipixno,irezid,istddev,iglat,iglon,iraw,isum,isel
|
|
parameter (TROWMAX=500000)
|
|
integer dpixno(TROWMAX),idx,sel(TROWMAX)
|
|
real dflux(TROWMAX),dglat(TROWMAX),dglon(TROWMAX)
|
|
real*8 dl,db
|
|
|
|
real dflux_49(TROWMAX),dflux_49_orig(TROWMAX),dflux_12(TROWMAX)
|
|
real dratio(TROWMAX)
|
|
integer*4 mask(TROWMAX)
|
|
|
|
! read DIRBE flux 4.9 microns (Resid)
|
|
! The sum of the selected column is 90891.360
|
|
! The mean of the selected column is 0.23114868
|
|
! The standard deviation of the selected column is 0.76384730
|
|
! The minimum of selected column is 3.64903696E-02
|
|
! The maximum of selected column is 110.34200
|
|
! The number of points used in calculation is 393216
|
|
!
|
|
! Expected value 0.0823
|
|
!
|
|
|
|
call ftgiou(unit,status)
|
|
call ftopen(unit,'/afs/mpa/project/integral/results/latscan/COBE/DIRBE_BAND04_ZSMA.FITS',0,block,status)
|
|
!call ftopen(unit,'/afs/mpa/project/integral/results/latscan/COBE/DIRBE_BAND3A_ZSMA.FITS',0,block,status)
|
|
if(status.ne.0) then
|
|
print*,'Bad file, skip it'
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
status=0
|
|
call exit(0)
|
|
endif
|
|
call FTGHDN(unit,nhdu)
|
|
call ftmahd(unit,2,hdutype,status)
|
|
call FTGHDN(unit,nhdu)
|
|
call FTGNRW(unit,tnrows, status)
|
|
call FTGNCL(unit,tncols, status)
|
|
|
|
if(tnrows.gt.trowmax) then
|
|
print*,'Too many points',tnrows,trowmax
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
stop
|
|
endif
|
|
exact=.false.
|
|
call ftgcno(unit,exact,'pixel_no',ipixno,status)
|
|
call ftgcno(unit,exact,'resid',irezid,status)
|
|
|
|
do iraw=1,tnrows
|
|
call ftgcvj(unit,ipixno,iraw,1,1,0,dpixno(iraw),anynull,status)
|
|
call ftgcve(unit,irezid,iraw,1,1,0.,dflux_49(iraw),anynull,status)
|
|
end do
|
|
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
if(status.gt.0) then
|
|
call exit(0)
|
|
endif
|
|
|
|
!
|
|
! Read DIRBE flux 1.25 microns
|
|
!
|
|
! The sum of the selected column is 194936.15
|
|
! The mean of the selected column is 0.49574827
|
|
! The standard deviation of the selected column is 1.3725306
|
|
! The minimum of selected column is 2.56427769E-02
|
|
! The maximum of selected column is 215.81805
|
|
! The number of points used in calculation is 393216
|
|
!
|
|
! Expected value: 0.0775
|
|
!
|
|
call ftgiou(unit,status)
|
|
call ftopen(unit,'/afs/mpa/project/integral/results/latscan/COBE/DIRBE_BAND1A_ZSMA.FITS',0,block,status)
|
|
if(status.ne.0) then
|
|
print*,'Bad file, skip it'
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
status=0
|
|
call exit(0)
|
|
endif
|
|
call FTGHDN(unit,nhdu)
|
|
call ftmahd(unit,2,hdutype,status)
|
|
call FTGHDN(unit,nhdu)
|
|
call FTGNRW(unit,tnrows, status)
|
|
call FTGNCL(unit,tncols, status)
|
|
|
|
if(tnrows.gt.trowmax) then
|
|
print*,'Too many points',tnrows,trowmax
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
stop
|
|
endif
|
|
exact=.false.
|
|
call ftgcno(unit,exact,'pixel_no',ipixno,status)
|
|
call ftgcno(unit,exact,'resid',irezid,status)
|
|
|
|
do iraw=1,tnrows
|
|
call ftgcvj(unit,ipixno,iraw,1,1,0,dpixno(iraw),anynull,status)
|
|
call ftgcve(unit,irezid,iraw,1,1,0.,dflux_12(iraw),anynull,status)
|
|
end do
|
|
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
if(status.gt.0) then
|
|
call exit(0)
|
|
endif
|
|
|
|
! read DIRBE coordinates
|
|
call ftgiou(unit,status)
|
|
call ftopen(unit,'/afs/mpa/project/integral/results/latscan/COBE/DIRBE_SKYMAP_INFO.FITS',0,block,status)
|
|
if(status.ne.0) then
|
|
print*,'Bad file, skip it'
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
status=0
|
|
call exit(0)
|
|
endif
|
|
call FTGHDN(unit,nhdu)
|
|
call ftmahd(unit,2,hdutype,status)
|
|
call FTGHDN(unit,nhdu)
|
|
call FTGNRW(unit,tnrows, status)
|
|
call FTGNCL(unit,tncols, status)
|
|
|
|
if(tnrows.gt.trowmax) then
|
|
print*,'Too many points',tnrows,trowmax
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
stop
|
|
endif
|
|
exact=.false.
|
|
call ftgcno(unit,exact,'QSPIXEL',ipixno,status)
|
|
call ftgcno(unit,exact,'GLON-CSC',iglon,status)
|
|
call ftgcno(unit,exact,'GLAT-CSC',iglat,status)
|
|
|
|
do iraw=1,tnrows
|
|
call ftgcvj(unit,ipixno,iraw,1,1,0,dpixno(iraw),anynull,status)
|
|
call ftgcve(unit,iglon,iraw,1,1,0.,dglon(iraw),anynull,status)
|
|
call ftgcve(unit,iglat,iraw,1,1,0.,dglat(iraw),anynull,status)
|
|
end do
|
|
|
|
call ftclos(unit,status)
|
|
call ftfiou(unit,status)
|
|
if(status.gt.0) then
|
|
call exit(0)
|
|
endif
|
|
|
|
! modify coordinates
|
|
do iraw=1,tnrows
|
|
if(dglon(iraw) .gt. 180.0) then
|
|
dglon(iraw)=dglon(iraw)-360.0
|
|
endif
|
|
enddo
|
|
|
|
! substruct background
|
|
do iraw=1,tnrows
|
|
dflux_12(iraw)=dflux_12(iraw)-0.07
|
|
dflux_49(iraw)=dflux_49(iraw)-0.08
|
|
enddo
|
|
|
|
! correction for absorbtion
|
|
do iraw=1,tnrows
|
|
dratio(iraw)=0.0
|
|
if(dflux_49(iraw).gt.0.5.and.dflux_49(iraw).lt.50.0.and.dflux_12(iraw).gt. 0.8 .and. dflux_12(iraw).lt.20.0) then
|
|
dratio(iraw)=dflux_12(iraw)/dflux_49(iraw)/3.5
|
|
else
|
|
dflux_49(iraw)=0.0
|
|
dflux_12(iraw)=0.0
|
|
dratio(iraw)=0.0
|
|
endif
|
|
enddo
|
|
|
|
do iraw=1,tnrows
|
|
if(dratio(iraw) .gt. 0.0 .and. dratio(iraw) .lt. 1.0) then
|
|
dratio(iraw)=(dratio(iraw))**(-0.25)
|
|
dflux_49(iraw)=dratio(iraw)*dflux_49(iraw)
|
|
else
|
|
dflux_49(iraw)=0.0
|
|
endif
|
|
enddo
|
|
|
|
do iraw=1,tnrows
|
|
if(abs(dglat(iraw)).gt.7.0) then
|
|
dflux_49(iraw)=0.0
|
|
endif
|
|
enddo
|
|
|
|
|
|
return
|
|
end
|
|
|