G.2 STACK_BIS.PRO

This program stacks a set of GAIA BBP simulated observations generated by SIM_STACK.PRO into a flux map using the alternative technique mentioned at the end of Section 6.4.

  1 ;   PROGRAM                 STACK_BIS.PRO
  2 
  3 ; Last Updated               07 Jun 2000           Mattia Vaccari
  4 
  5 ; This program stacks a set of GAIA BBP simulated observations of an HST WFPC2
  6 ;  field generated by SIM_STACK.PRO into a GAIA flux map.
  7 
  8 lsc:
  9 print,'Do you want to display GAIA flux map on screen (y/n)?'
 10 sc='...' & read,sc
 11 if sc ne 'y' and sc ne 'n' then begin
 12     print,'Invalid choice!'
 13     goto,lsc
 14     endif
 15 
 16 lps:
 17 print,'Do you want to save the flux map in magnitudes to an eps file (y/n)?'
 18 ps='...' & read,ps
 19 if ps ne 'y' and ps ne 'n' then begin
 20     print,'Invalid choice!'
 21     goto,lps
 22     endif
 23 if ps eq 'y' then begin
 24     print,'Enter filename without the .eps extension:'
 25     fn='fluxmap.eps' & read,fn
 26     print,'Enter the side of the image in centimetres:'
 27     fms=16. & read,fms
 28     endif
 29 
 30 lst:
 31 print,'Do you want to save the important data into an IDL structure (y/n)?'
 32 st='...' & read,st
 33 if st ne 'y' and st ne 'n' then begin
 34     print,'Invalid choice!'
 35     goto,lst
 36     endif
 37 if st eq 'y' then begin
 38     lstfile:
 39     print,'Do you want to save this structure into a dat file (y/n)?'
 40     stfile='...' & read,stfile
 41     if stfile ne 'y' and stfile ne 'n' then begin
 42         print,'Invalid choice!'
 43         goto,lstfile
 44         endif
 45     if stfile eq 'y' then begin
 46         print,'Enter filename without the .dat extension:'
 47         simfn='sim_data.dat' & read,simfn
 48         endif
 49     endif
 50 
 51 ; RETRIEVE DATA FROM SIMSTR STRUCTURE
 52 dadsfile=simstr.dadsfile & date=simstr.date & instrume=simstr.instrume
 53 camerastr=simstr.camerastr &targname=simstr.targname & ra_targ=simstr.ra_targ
 54 dec_targ=simstr.dec_targ & pa_v3=simstr.pa_v3 & filtnam1=simstr.filtnam1
 55 filtnam2=simstr.filtnam2 & filtnu=simstr.filtnu & atodgain=simstr.atodgain
 56 gain=simstr.gain & photflam=simstr.photflam & photzpt=simstr.photzpt
 57 calpar=simstr.calpar & vi=simstr.vi & dm=simstr.dm & hpix=simstr.hpix
 58 exptime=simstr.exptime & hstdnc=simstr.hstdnc
 59 ; hstelc=simstr.hstelc & hstmagins=simstr.hstmagins
 60 muback=simstr.muback & parhstdn=simstr.parhstdn & hstmagsta=simstr.hstmagsta
 61 hstmag=simstr.hstmag & rdim=simstr.rdim & hstmagcen=simstr.hstmagcen
 62 ss=simstr.ss & ssstr=simstr.ssstr & nobs=simstr.nobs & sd=simstr.sd
 63 posa=simstr.posa & cx=simstr.cx & cy=simstr.cy & gaiaet=simstr.gaiaet
 64 rn=simstr.rn & spfx=simstr.spfx & spfy=simstr.spfy & obs=simstr.obs
 65 ssfx=simstr.ssfx & ssfy=simstr.ssfy & gstep=simstr.gstep & gsize=simstr.gsize
 66 parn=simstr.parn & px=simstr.px & py=simstr.py & psx=simstr.psx
 67 psy=simstr.psy & sx=simstr.sx & sy=simstr.sy & xi=simstr.xi & yi=simstr.yi
 68 
 69 print,'Stacking GAIA ',ssstr,' pixels/sample BBP simulated observations.'
 70 print,'Please wait...'
 71 
 72 ; FLUX MAP GRID IN THE HST FOR
 73 fmgrid=fltarr(2,gsize^2)
 74 k=0l
 75 for i=0,gsize-1 do for j=0,gsize-1 do begin
 76     fmgrid(*,k)=[i-(gsize-1)/2,j-(gsize-1)/2]*gstep
 77     k=k+1
 78     endfor
 79 ; fmgrid is a (2,gsize^2) array containing along its 2 columns the x and y
 80 ;  coordinates of the gsize^2 flux map grid points in the HST FOR.
 81 ;  The first row contains the coordinate of the point at the lower left corner
 82 ;  of the grid and the two columns are ordered first bottom/up and then
 83 ;  left/right (which is the order in which IDL dipslays 2-D arrays).
 84 ;  The essentially 2-D arrays folded into 1-D arrays that will be used in the
 85 ;  following for the sake of efficiency will maintain the same ordering.
 86 
 87 ; FLUX MAP GRID IN THE OBSERVATION FOR
 88 ones=replicate(1.,gsize^2)
 89 ;obsfmx=[[+cos(posa)],[+sin(posa)]]#fmgrid+$
 90 ; diag([[+cos(posa)],[+sin(posa)]]#transpose([[cx],[cy]]))#ones
 91 ;obsfmy=[[-sin(posa)],[+cos(posa)]]#fmgrid+$
 92 ; diag([[-sin(posa)],[+cos(posa)]]#transpose([[cx],[cy]]))#ones
 93 ;obsfmx=[[+cos(posa)],[-sin(posa)]]#fmgrid+cx#transpose(ones)
 94 ;obsfmy=[[+sin(posa)],[+cos(posa)]]#fmgrid+cy#transpose(ones)
 95 obsfmx=[[+cos(posa)],[+sin(posa)]]#fmgrid-$
 96 diag([[+cos(posa)],[+sin(posa)]]#transpose([[cx],[cy]]))#ones
 97 obsfmy=[[-sin(posa)],[+cos(posa)]]#fmgrid-$
 98 diag([[-sin(posa)],[+cos(posa)]]#transpose([[cx],[cy]]))#ones
 99 ; obsfmx and obsfmy are (nobs,gsize^2) arrays containing the x and y
100 ;  coordinates of the flux map elements in the observation for
101 
102 ; OBSERVATIONS' REBINNING
103 obsreb=rebin(obs,nobs,xi*ssfx,yi*ssfy,/sample)/(ssfx*ssfy)
104 
105 ; ASSIGNMENT OF A SUBSAMPLE VALUE TO EACH FLUX MAP ELEMENT
106 xin=round(obsfmx/(sx/ssfx)+0.5*(xi*ssfx-1))
107 yin=round(obsfmy/(sy/ssfy)+0.5*(yi*ssfy-1))
108 ; xin and yin express obsfmx and obsfmy in samples
109 nfm=fltarr(nobs,gsize^2)
110 for n=0,nobs-1 do nfm(n,*)=obsreb([intarr(gsize^2)+n],[xin(n,*)],[yin(n,*)])
111 ;nfm=obsreb(indgen(nobs)#replicate(1,gsize^2),$
112 ; round(obsfmx/(sx/ssfx)+0.5*(xi*ssfx-1)),$
113 ; round(obsfmy/(sy/ssfy)+0.5*(yi*ssfy-1)))
114 ; nfm now contains along its nobs columns the sample value of the sample
115 ;  whose centre is nearest to each flux map element for each observation
116 ; In other words, nfm contains along its nobs columns the nobs
117 ;  one-dimensional partial flux maps
118 
119 nfm=transpose(reform(nfm,nobs,gsize,gsize),[0,2,1])
120 ; nfm now contains the nobs two-dimensional partial flux maps
121 gfm=total(nfm,1)
122 ; gfm is the two-dimensional GAIA glaobal flux map
123 
124 ; FLUX MAP IN MAGNITUDES
125 gfmmag=-2.5*alog10(gfm>parn)
126 gfmmax=max(gfmmag,min=gfmmin)
127 ; GAIA flux map in uncalibrated magnitudes
128 pargfmdn=parn/gain
129 gfmdnc=gfm/gain
130 ; GAIA flux map in data number counts
131 gfmmagsta=-2.5*alog10((gfmdnc>pargfmdn)/(gaiaet*nobs*(sx*sy/1000.)^2))+$
132 calpar(0,filtnu)+calpar(1,filtnu)*vi+calpar(2,filtnu)*vi^2+dm+0.1
133 ; GAIA flux map in stamag standard magnitudes
134 
135 ; FLUX MAP ON SCREEN
136 ;if sc eq 'y' then begin
137 ;mag=2.
138 ;dispsize=mag*gsize
139 ;set_plot,'x'
140 ;loadct,0
141 ;window,!d.window+1,xsize=dispsize,ysize=dispsize
142 ;pfx=[0,1,1,0] & pfy=[0,0,1,1]
143 ;for i=0,gsize-1 do for j=0,gsize-1 do begin
144 ;polyfill,mag*(pfx+i),mag*(pfy+j),color=round((gfmmag(i,j)-gfmmin)/$
145 ; (gfmmax-gfmmin)*255.),/device
146 ;endfor
147 ;endif
148 
149 if sc eq 'y' then begin
150     mag=2.
151     dispsize=mag*gsize
152     gfmmagsc=rebin(gfmmag,dispsize,dispsize,/sample)
153     set_plot,'x'
154     loadct,0
155     window,!d.window+1,xsize=dispsize,ysize=dispsize
156     tvscl,gfmmagsc
157     endif
158 
159 ; SAVE FLUX MAP INTO AN EPS FILE
160 if ps eq 'y' then begin
161     set_plot,'ps'
162     device,filename=fn+'.eps',xsize=fms,ysize=fms,xoffset=(21.-fms)/2.,$
163     yoffset=(29.7-fms)/2.,/encapsulated,bits_per_pixel=8
164     tvscl,gfmmag,xsize=fms,ysize=fms,/centimeters
165     device,/close
166     set_plot,'x'
167     endif
168 ; Note that the image is centered on an a4 page
169 
170 ; SAVE IMPORTANT DATA INTO AN IDL STRUCTURE
171 if st eq 'y' then begin
172     ; Simulation Structure
173     result=execute('simstr_bis={dadsfile:dadsfile,date:date,instrume:instrume'+$
174     ',camerastr:camerastr,targname:targname,ra_targ:ra_targ,dec_targ:dec_targ'+$
175     ',pa_v3:pa_v3,filtnam1:filtnam1,filtnam2:filtnam2,filtnu:filtnu'+$
176     ',atodgain:atodgain,gain:gain,photflam:photflam,photzpt:photzpt'+$
177     ',calpar:calpar,vi:vi,dm:dm,hpix:hpix,exptime:exptime,hstdnc:hstdnc'+$
178     ; ',hstelc:hstelc,hstmagins:hstmagins'+$
179     ',muback:muback,parhstdn:parhstdn,hstmagsta:hstmagsta,hstmag:hstmag'+$
180     ',rdim:rdim,hstmagcen:hstmagcen'+$
181     ',ss:ss,ssstr:ssstr,nobs:nobs,sd:sd,posa:posa,cx:cx,cy:cy,gaiaet:gaiaet'+$
182     ',rn:rn,spfx:spfx,spfy:spfy,obs:obs'+$
183     ',ssfx:ssfx,ssfy:ssfy,gstep:gstep,gsize:gsize,nfm:nfm,gfm:gfm,parn:parn'+$
184     ',gfmmag:gfmmag,gfmmagsta:gfmmagsta'+$
185     ',px:px,py:py,psx:psx,psy:psy,sx:sx,sy:sy,xi:xi,yi:yi}')
186     endif
187 
188 ; SAVE IMPORTANT DATA INTO A DAT FILE
189 if stfile eq 'y' then begin
190     save,filename=simfn+'.dat',/verbose,simstr_bis
191     endif
192 
193 end