program cal17_1 *********************************************************************** * create calibration files (336) * * January 30, 2007 H. Koshiishi * *********************************************************************** c integer naxes(2),nantst(84),icimg(1000),iimg(1000), - pir(2000),pil(2000),plr(2000),pll(2000),pwr(2000),pwl(2000) integer*2 iccald(336) real gbeam(512,512,0:10),flxr(84),flxl(84), - mapra(512,512),maprb(512,512),mapla(512,512),maplb(512,512), - cmapr(512,512),cmapl(512,512),pmat(4), - par(2000),pal(2000) real rlph(84) character*80 filehd,infpar,outf,outfft character*18 dattyp character*12 dateframe character*1 fform c ccc write(6,'(/'' input directory for output files : ''$)') read(5,'(a)') outf ccc write(6,'('' input rawdata filename header : ''$)') read(5,'(a)') filehd ccc write(6,'('' input parameter filename : ''$)') read(5,'(a)') infpar c ccc write(6,'('' input first frame number : ''$)') read(5,*) ncfrst if (ncfrst.ne.0) then ccc write(6,'('' input last frame number : ''$)') read(5,*) ncfrend ccc write(6,'('' input frame interval : ''$)') read(5,*) ictvl ncimg=0 do 150 i=ncfrst,ncfrend,ictvl ncimg=ncimg+1 icimg(ncimg)=i 150 continue else ccc write(6,'('' input number of frames : ''$)') read(5,*) ncimg ccc write(6,'('' input frame numbers : ''$)') read(5,*) (icimg(i),i=1,ncimg) end if ccc write(6,'('' input integration frames'')') ccc write(6,'('' for calibration : ''$)') read(5,*) iccal ccc write(6,'('' for image processing : ''$)') read(5,*) ictgr c ccc write(6,'('' quiet : 1 / burst : 2 : ''$)') read(5,*) iqb c ccc write(6,'('' input pixel-size of the output files : ''$)') read(5,*) npix c ccc write(6,'('' input offset of the image center'')') ccc write(6,'('' unit(1:arcsec/2:pixel) & x-/y-offsets : ''$)') read(5,*) nuoff,xoffa,yoffa if (nuoff .eq. 1) then xoffp=xoffa/4.91106 yoffp=yoffa/4.91106 else xoffp=xoffa yoffp=yoffa xoffa=xoffp*4.91106 yoffa=yoffp*4.91106 endif c ccc write(6,'('' input first frame number : ''$)') read(5,*) nfrst if (nfrst.ne.0) then ccc write(6,'('' input last frame number : ''$)') read(5,*) nfrend ccc write(6,'('' input frame interval : ''$)') read(5,*) itvl nimg=0 do 160 i=nfrst,nfrend,itvl nimg=nimg+1 iimg(nimg)=i 160 continue else ccc write(6,'('' input number of frames : ''$)') read(5,*) nimg ccc write(6,'('' input frame numbers : ''$)') read(5,*) (iimg(i),i=1,nimg) end if ccc write(6,'('' input integration frames'')') ccc write(6,'('' for calibration : ''$)') read(5,*) ical ccc write(6,'('' for image processing : ''$)') read(5,*) itgr c ccc write(6,'('' input number limit of criterion(r+l,r-l): ''$)') read(5,*) nmaxa,nmaxs ccc write(6,'('' input lower limit of criterion(r+l,r-l): ''$)') read(5,*) cradd,crsub ccc write(6,'('' input factor of criterion : ''$)') read(5,*) crfac c ***** default parameters ***** c loutf=index(outf,' ')-1 c call cft512(mapra,maprb,0) c iunit=11 c call rdpar(iunit,infpar, - dec1,dec2,dec3,ha1,ha2,ha3, - solr1,solr2,solr3,solp1,solp2,solp3,solb1,solb2,solb3, - nantst) c call msrc(nantst,1.,gbeam) c iounit=2 nbit=32 naxis=2 naxes(1)=336 naxes(2)=1 c do 100 iifr=1,ncimg c iframe=icimg(iifr) c ***** initialize rlph ***** c if (iifr.eq.1) rlph(1)=1. c call snapuv17(iunit,filehd,iframe,iframe+ictgr-1, - nantst,iccal,iccald,rlph,1., - ndjst,msjst,msjsts,msjste,natt,nfalt,nfreq,nstat, - flxr,flxl,mapra,maprb,mapla,maplb,nfrmst) if (nfrmst.ne.0) go to 100 c call cft512(mapra,maprb,-1) call cft512(mapla,maplb,-1) c call caleph(dec1,dec2,dec3,ha1,ha2,ha3, - solr1,solr2,solr3,solp1,solp2,solp3,solb1,solb2,solb3, - msjst,dec,ha,solr,solp,solb,az,alt,za,pmat) c ***** preclean process for brightness normalization ***** c cradd=3000. crsub=3000. crfac=0. nmaxa=2000 nmaxs=2000 c do 200 l=1,512 do 200 i=1,512 200 maprb(i,l)=mapra(i,l) c call clean17(3,maprb,pmat,solr*1.0125, - nantst,gbeam,cradd,crfac,cra,0.,0.,corf, - nmaxa,bgrr,dskr,pir,plr,pwr,par,cmapr) c do 210 l=1,512 do 210 i=1,512 210 maprb(i,l)=mapla(i,l) c call clean17(3,maprb,pmat,solr*1.0125, - nantst,gbeam,cradd,crfac,cra,0.,0.,corf, - nmaxs,bgrl,dskl,pil,pll,pwl,pal,cmapl) c ***** add and subtract rcp/lcp ***** c do 300 l=1,512 do 300 i=1,512 300 maprb(i,l)=mapra(i,l)-bgrr-(mapla(i,l)-bgrl)/dskl*dskr do 310 l=1,512 do 310 i=1,512 310 mapra(i,l)=(mapra(i,l)+mapla(i,l)/dskl*dskr)/2. do 320 l=1,512 do 320 i=1,512 320 mapla(i,l)=maprb(i,l)/2. c ***** clean process ***** c pxew=999. pxns=999. nmaxa=2000 nmaxs=2000 c call clean17(2,mapra,pmat,solr*1.0125, - nantst,gbeam,cradd,crfac,cra,pxew,pxns,corf, - nmaxa,bgr,dskbr,pir,plr,pwr,par,cmapr) c lfilehd=index(filehd,' ')-1 if (iframe.lt.10) - write(dateframe,'(a6,''_0000'',i1)') - filehd(lfilehd-5:lfilehd),iframe if (iframe.ge.10 .and. iframe.lt.100) - write(dateframe,'(a6,''_000'',i2)') - filehd(lfilehd-5:lfilehd),iframe if (iframe.ge.100 .and. iframe.lt.1000) - write(dateframe,'(a6,''_00'',i3)') - filehd(lfilehd-5:lfilehd),iframe if (iframe.ge.1000 .and. iframe.lt.10000) - write(dateframe,'(a6,''_0'',i4)') - filehd(lfilehd-5:lfilehd),iframe if (iframe.ge.10000) - write(dateframe,'(a6,''_'',i5)') - filehd(lfilehd-5:lfilehd),iframe c fform='f' c outfft=outf(1:loutf)//'/i'//fform//'c'//dateframe call ftinit(iounit,outfft,2880,istat) if (istat.ne.0) go to 900 dattyp='calibration_data' call puthdr(iounit,nbit,naxis,naxes, - ndjst,msjst,msjsts,msjste,iframe,iframe+ictgr-1, - 0,natt,nfalt,nfreq,nstat,dattyp) c call putcal(iounit,iccal,ictgr,pxew,pxns,corf) c call ftpdef(iounit,nbit,naxis,naxes,0,0,istat) call ftppri(iounit,0,1,naxes(1)*naxes(2),iccald,istat) call ftclos(iounit,istat) c if (istat.ne.0) go to 900 write(6,'('' out file :'',a40)') outfft c write(6,*) 'pxew , pxns =',pxew,',',pxns c 100 continue c 900 write(6,'('' istat ='',i8)') istat c end c c *********************************************************************** * putcal - write calibration parameters to fits file * * iunit : input,integer fortran unit number for fits file * * nfrcal : input,integer number of calibration frames * * nfrtgr : input,integer number of integration frames * * pxew : input,real x-offset of the center of dirty disk * * pxns : input,real y-offset of the center of dirty disk * * corf : input,real correlation of dirty disk and model * * January 30, 2007 H. Koshiishi * *********************************************************************** c subroutine putcal(iunit,nfrcal,nfrtgr,pxew,pxns,corf) c integer nfrcal,nfrtgr real pxew,pxns,corf character*44 cmnt c cmnt=' ' c call ftpkys(iunit,'progname','cal17_1 H.Koshiishi',cmnt,istat) call ftpkyj(iunit,'nfrcal',nfrcal,'number of calibration frames', - istat) call ftpkyj(iunit,'nfrtgr',nfrtgr,'number of integration frames', - istat) call ftpkye(iunit,'ddoff1',pxew,8, - 'x-offset of the dirty disk',istat) call ftpkye(iunit,'ddoff2',pxns,8, - 'y-offset of the dirty disk',istat) call ftpkyf(iunit,'ddcorr',corf,4, - 'correlation between dirty disk and model',istat) c if (istat.ne.0) write(6,'('' putcal : istat ='',i6)') istat return end