| 1 | #####################################################################################
|
|---|
| 2 | #### Commands to run the different programs to produce foreground maps
|
|---|
| 3 | #### and compute radio-source subtracted P(k)
|
|---|
| 4 | #####################################################################################
|
|---|
| 5 |
|
|---|
| 6 | ### Cube definition in file cubedef.h
|
|---|
| 7 |
|
|---|
| 8 | ### Step 1/ Produce an LSS data cube with appropriate size and redshift using SimLSS
|
|---|
| 9 | # 1.a/ Run SimLSS
|
|---|
| 10 | csh> ~/Objs/exe/cmvginit3df -a -1 -2 -C -G 0. -F 0 -x 360,3 -y 360,3 -z 256,1.5 -Z 0.56 -8 1. -n 10000 -O 0,2 -o lssz056 -T 2
|
|---|
| 11 | # 1.b/ To run SimLSS with GSM map parameters (DeltaFreq=500 MHz)
|
|---|
| 12 | csh> ~/Objs/exe/cmvginit3df -a -1 -2 -C -G 0. -F 0 -x 360,3 -y 360,3 -z 256,3 -Z 0.60 -8 1. -n 10000 -O 0,2 -o lssz060 -T 2
|
|---|
| 13 |
|
|---|
| 14 | # 1.c/ To run SimLSS with GSM map parametersand 40 (phi/alpha) x 30 (theta/delta) deg maps (DeltaFreq=500 MHz) @ z=0.6
|
|---|
| 15 | csh> ~/Objs/exe/cmvginit3df -a -1 -2 -C -G 0. -F 0 -x 600,1.9 -y 800,1.9 -z 256,2.8 -Z 0.60 -8 1. -n 10000 -O 0,2 -o lssz060 -T 2
|
|---|
| 16 |
|
|---|
| 17 | # 1.d/ To run SimLSS with GSM map parametersand 90x30 deg maps (DeltaFreq=500 MHz) @ z=1 [ 90deg-> phi/alpha, 30deg -> theta/delta]
|
|---|
| 18 | csh> ~/Objs/exe/cmvginit3df -a -1 -2 -C -G 0. -F 0 -x 600,2.9 -y 1800,2.9 -z 256,3.5 -Z 1.0 -8 1. -n 10000 -O 0,2 -o lssz100 -T 2
|
|---|
| 19 |
|
|---|
| 20 | csh> ~/Objs/exe/cmvginit3df -a -1 -2 -C -G 0. -F 0 -x 600,3.8 -y 1800,3.8 -z 256,4.2 -Z 1.5 -8 1. -n 10000 -O 0,2 -o lssz150 -T 2
|
|---|
| 21 |
|
|---|
| 22 | # 1.c/ Change the X and Z axis of the cube to adapt it to RadioBeam package convention
|
|---|
| 23 | # SimLSS output : the radial (redshift) direction along X axis of the cube (TArray)
|
|---|
| 24 | # RadioBeam cubes : the radial (redshift) direction along Z axis of the cube (TArray)
|
|---|
| 25 | # Execucte the following script in spiapp :
|
|---|
| 26 |
|
|---|
| 27 | csh> cat > racube.pic
|
|---|
| 28 | set f lssz060
|
|---|
| 29 | readfits ${f}_r.fits
|
|---|
| 30 | rename ${f}_r map
|
|---|
| 31 | print map
|
|---|
| 32 | c++exec \
|
|---|
| 33 | TArray<r_4> omap(map.SizeY(),map.SizeZ(),map.SizeX()-2 ); \
|
|---|
| 34 | for(sa_size_t i=0;i<omap.SizeX();i++) \
|
|---|
| 35 | for(sa_size_t j=0;j<omap.SizeY();j++) \
|
|---|
| 36 | for(sa_size_t k=0;k<omap.SizeZ();k++) \
|
|---|
| 37 | omap(i,j,k)=map(k+1,i,j); \
|
|---|
| 38 | KeepObj(omap);
|
|---|
| 39 |
|
|---|
| 40 | rename omap lsscube
|
|---|
| 41 | print lsscube
|
|---|
| 42 | # expmeansig lsscube val
|
|---|
| 43 | saveppf lsscube lsscubez060.ppf
|
|---|
| 44 |
|
|---|
| 45 | csh> spiapp -term -exec racube.pic
|
|---|
| 46 |
|
|---|
| 47 | #### Cube LSS 40x30 deg (3') @ z=0.6 ( lsscubez060.ppf )
|
|---|
| 48 | #### -> Size= 122880000 Mean=-7.01664e-05 Sigma=2.53016 Min=-13.7439 Max=14.4648
|
|---|
| 49 |
|
|---|
| 50 | ## Step 2/ Produce synchrotron and radio source sky cubes (cube unit is Temparature- Kelvin)
|
|---|
| 51 | # 2.a/ Synchrotron map from HASLAM 400 MHz map
|
|---|
| 52 | csh> ./Objs/syncube syncmap_eq.fits syncube.ppf syncmap.ppf
|
|---|
| 53 | # 2.b/ radio source cube from NVSS catalog
|
|---|
| 54 | csh> ./Objs/srcat2cube -nvss nvss.fits nvsscube.ppf nvssmap.ppf
|
|---|
| 55 | # Or from the north20 catalog :
|
|---|
| 56 | csh> ./Objs/srcat2cube -north20 north20cm.fits north20cube.ppf north20map.ppf
|
|---|
| 57 |
|
|---|
| 58 | # 2.c/ Add the two cubes using the following spiapp script
|
|---|
| 59 | csh> cat > sumcubes.pic
|
|---|
| 60 | openppf syncube.ppf
|
|---|
| 61 | openppf radsrccube.ppf
|
|---|
| 62 | # expmeansig syncube val
|
|---|
| 63 | # expmeansig nvsscube val
|
|---|
| 64 | c++exec TArray<r_4> fgndcube = syncube+radsrccube; KeepObj(fgndcube);
|
|---|
| 65 | print fgndcube
|
|---|
| 66 | # expmeansig fgndcube val
|
|---|
| 67 | saveppf fgndcube fgndcube.ppf
|
|---|
| 68 |
|
|---|
| 69 | csh> spiapp -term -exec sumcubes.pic
|
|---|
| 70 |
|
|---|
| 71 | #### syncube:Mean= 1.8101 Sigma= 0.326538 Min= 0.857019 Max= 3.58987
|
|---|
| 72 | #### nvsscube: Mean= 1.95073 Sigma= 1.68515 Min= 0.857019 Max= 428.398
|
|---|
| 73 | #### fgndcube=syncube+nvsscube: Mean= 0.140623 Sigma= 1.65068 Min= 0 Max= 426.559
|
|---|
| 74 | #### north20: fgndcube_north:
|
|---|
| 75 |
|
|---|
| 76 | ## Step 2.b/ Produce foreground cube from GSM
|
|---|
| 77 | csh> ./Objs/gsm2cube ../Catalogs/GSM/ 1 256 fgndcube_gsm.ppf
|
|---|
| 78 |
|
|---|
| 79 | ## Step 3/ Apply lobe (50 meter diameter array) effect on foreground cube and LSS cube
|
|---|
| 80 | csh> set ddish=55.
|
|---|
| 81 | csh> set ddishcor=55.
|
|---|
| 82 | csh> ./Objs/applobe $ddish fgndcube.ppf fgndcube_lobe.ppf
|
|---|
| 83 | csh> ./Objs/applobe -fib $ddish fgndcube.ppf fgndcube_flobe.ppf
|
|---|
| 84 | csh> ./Objs/applobe $ddish lsscube.ppf lsscube_lobe.ppf
|
|---|
| 85 | csh> ./Objs/applobe -fib $ddish lsscube.ppf lsscube_flobe.ppf
|
|---|
| 86 | ## Step 3.b/ Correct for the lobe effect by bringing all to the beam of Diam/Lambda = 150 (55 m @ z=0.7 - 820 MHz)
|
|---|
| 87 | csh> ./Objs/applobe $ddish lsscube_lobe.ppf lsscube_corlobe.ppf $ddishcor
|
|---|
| 88 | csh> ./Objs/applobe $ddish fgndcube_lobe.ppf fgndcube_corlobe.ppf $ddishcor
|
|---|
| 89 |
|
|---|
| 90 | ## Step 3.c/ Apply lobe (Filled 11x11 5m dishes array) effect on foreground cube and LSS cube
|
|---|
| 91 | csh> ./Objs/applobe repf11x11.ppf fgndcube.ppf fgndcube_lobe.ppf
|
|---|
| 92 | csh> ./Objs/applobe repf11x11.ppf lsscube.ppf lsscube_lobe.ppf
|
|---|
| 93 | ## Step 3.d/ Correct for the lobe effect by bringing all to the beam of Diam/Lambda = 150 (55 m @ z=0.7 - 820 MHz)
|
|---|
| 94 | csh> ./Objs/applobe repf11x11.ppf lsscube_lobe.ppf lsscube_corlobe.ppf $ddishcor
|
|---|
| 95 | csh> ./Objs/applobe repf11x11.ppf fgndcube_lobe.ppf fgndcube_corlobe.ppf $ddishcor
|
|---|
| 96 |
|
|---|
| 97 | ### Step 4/ Compute power spectra
|
|---|
| 98 | ## mass to temperature converion factor CT21 ~= 0.21 mK for gHI=2% , 0.11 for gHI=1% , 0.13 for gHI=0.008x(1+0.6)
|
|---|
| 99 | ## Foreground maps are in temperature
|
|---|
| 100 | ## Noise fluctuations Sigma^2 ~ T_sys^2 / t_obs * DeltaFreq
|
|---|
| 101 | ## Tsys ~ 50 K , DeltaFreq ~ 0.5 MHz , t_obs ~ 1 day ~ 80 000 s.
|
|---|
| 102 | ## sigma_noise ~ 0.25 mK -> 3 mK
|
|---|
| 103 | # 4.a/ LSS power spectrum without noise
|
|---|
| 104 | csh> ./Objs/calcpk lsscube.ppf lsspk.ppf 0.13
|
|---|
| 105 | # and with noise
|
|---|
| 106 | csh> ./Objs/calcpk lsscube.ppf lsspkwn.ppf 0.13 3
|
|---|
| 107 | # with the lobe effect
|
|---|
| 108 | csh> ./Objs/calcpk lsscube_lobe.ppf lsspklobe.ppf 0.13
|
|---|
| 109 | csh> ./Objs/calcpk lsscube_flobe.ppf lsspkflobe.ppf 0.13
|
|---|
| 110 | csh> ./Objs/calcpk lsscube_lobe.ppf lsspklobewn.ppf 0.13 3
|
|---|
| 111 | csh> ./Objs/calcpk lsscube_corlobe.ppf lsspkcorlobe.ppf 0.13
|
|---|
| 112 |
|
|---|
| 113 | # 4.b/ Foreground power spectrum
|
|---|
| 114 | csh> ./Objs/calcpk fgndcube.ppf fgndpk.ppf 1000
|
|---|
| 115 | csh> ./Objs/calcpk fgndcube_lobe.ppf fgndpklobe.ppf 1000
|
|---|
| 116 | csh> ./Objs/calcpk fgndcube_flobe.ppf fgndpkflobe.ppf 1000
|
|---|
| 117 | csh> ./Objs/calcpk fgndcube_corlobe.ppf fgndpkcorlobe.ppf 1000
|
|---|
| 118 |
|
|---|
| 119 | # 4.c/ Extract LSS P(k) from Foreground+LSS+noise , after cleaning/subtraction without beam
|
|---|
| 120 | csh> set beamdesc=repf11x11.ppf
|
|---|
| 121 | csh> set ddishcor=55.
|
|---|
| 122 | csh> set noiselev=1.
|
|---|
| 123 | csh> ./Objs/calcpk2 lsscube.ppf 0.13 fgndcube.ppf 1000 subpk.ppf $noiselev $beamdesc 0. 0. P2
|
|---|
| 124 | # 4.d / Extract LSS P(k) from Foreground+LSS+noise and beam effect, without beam correction
|
|---|
| 125 | csh> ./Objs/calcpk2 lsscube_lobe.ppf 0.13 fgndcube_lobe.ppf 1000 subpklobe.ppf $noiselev $beamdesc 0. 0. P2
|
|---|
| 126 | # 4.e / Extract LSS P(k) from Foreground+LSS+noise and beam effect - correcting to a beam of Diam= $ddishcor
|
|---|
| 127 | csh> ./Objs/calcpk2 lsscube_lobe.ppf 0.13 fgndcube_lobe.ppf 1000 subpkcorlobe.ppf $noiselev $beamdesc $ddishcor 0. P2 reclsscorlobe.ppf
|
|---|
| 128 | # Or using a linear fit for foreground subtraction (old version)
|
|---|
| 129 | csh> ./Objs/calcpk2 lsscube_lobe.ppf 0.13 fgndcube_lobe.ppf 1000 subpkcorlobep1.ppf $noiselev $beamdesc $ddishcor 0. P2
|
|---|
| 130 | # 4.f / Estimate residual noise from Foreground removal :
|
|---|
| 131 | csh> ./Objs/calcpk2 lsscube.ppf 0. fgndcube_lobe.ppf 1000 residcorlobe.ppf $noiselev $beamdesc $ddishcor 0. P2
|
|---|
| 132 | csh> ./Objs/calcpk2 lsscube.ppf 0. fgndcube_lobe.ppf 1000 residnocor.ppf $noiselev $beamdesc 0. 0. P2
|
|---|
| 133 |
|
|---|
| 134 | ### Step 5 / Check the results using spiapp
|
|---|
| 135 | setaxesatt 'font=helvetica,bold,16 fixedfontsize minorticks'
|
|---|
| 136 | delobjs *
|
|---|
| 137 | openppf fgndpk.ppf
|
|---|
| 138 | openppf fgndpklobe.ppf
|
|---|
| 139 | openppf fgndpkflobe.ppf
|
|---|
| 140 | openppf fgndpkcorlobe.ppf
|
|---|
| 141 | openppf lsspk.ppf
|
|---|
| 142 | openppf lsspklobe.ppf
|
|---|
| 143 | openppf lsspkcorlobe.ppf
|
|---|
| 144 |
|
|---|
| 145 | openppf lsspkflobe.ppf
|
|---|
| 146 | openppf lsspklobewn.ppf
|
|---|
| 147 | openppf subpklobe.ppf
|
|---|
| 148 |
|
|---|
| 149 | openppf subpkcorlobe.ppf
|
|---|
| 150 |
|
|---|
| 151 | openppf residcorlobe.ppf
|
|---|
| 152 | openppf residnocor.ppf
|
|---|
| 153 | # openppf subpknolssnocor.ppf
|
|---|
| 154 |
|
|---|
| 155 | disp lsspk 'logx logy nsta xylimits=0.01,2.,4e-11,8e-6 gold'
|
|---|
| 156 | disp lsspklobe 'same nsta orange'
|
|---|
| 157 | disp lsspklobewn 'same nsta siennared'
|
|---|
| 158 | settitle ' Pk[LSS] - without normalisation' ' ' 'font=helvetica,bold,16 black'
|
|---|
| 159 |
|
|---|
| 160 |
|
|---|
| 161 | disp fgndpk 'logx logy nsta xylimits=0.01,2.,1e-10,1. navyblue'
|
|---|
| 162 | disp fgndpklobe 'same nsta blue'
|
|---|
| 163 | disp fgndpkcorlobe 'same nsta skyblue'
|
|---|
| 164 | disp lsspk 'same nsta gold'
|
|---|
| 165 | disp lsspkflobe 'same nsta yellow'
|
|---|
| 166 | disp subpkcorlobe 'same nsta red'
|
|---|
| 167 | disp residcorlobe 'same nsta green'
|
|---|
| 168 | disp residnocor 'same nsta forestgreen'
|
|---|
| 169 |
|
|---|
| 170 | disp lsspklobewn 'same nsta siennared'
|
|---|
| 171 | # settitle 'Pk[LSS] , Pk[Foreground] and lobe effect (Dish D=50 m)' ' ' 'font=helvetica,bold,18'
|
|---|
| 172 | settitle 'Pk[LSS] , Pk[Foreground=GSM] and lobe effect (Dish D=50 m)' ' ' 'font=helvetica,bold,18'
|
|---|
| 173 |
|
|---|
| 174 | set lines ( 'Pk[Foreground]' 'Pk[fgnd]*Lobe' 'Pk[fgnd]*Lobe/Corrected' 'Pk[LSS]' 'Pk[LSS]*Lobe+Noise' )
|
|---|
| 175 | set cols ( navyblue blue skyblue gold siennared )
|
|---|
| 176 | textdrawer lines cols 'font=helvetica,bold,16 frame'
|
|---|
| 177 |
|
|---|
| 178 |
|
|---|
| 179 | disp lsspk 'logx logy nsta xylimits=0.005,2.,4e-9,4e-5 gold'
|
|---|
| 180 | disp lsspklobewn 'same nsta siennared'
|
|---|
| 181 | disp subpkcorlobe 'same nsta red'
|
|---|
| 182 | disp subpknolss 'same nsta green'
|
|---|
| 183 |
|
|---|
| 184 | # Calcul du volume total en Mpc^3
|
|---|
| 185 | set VOL (1.9*1.9*2.8*800*600*256)
|
|---|
| 186 | # set VOL (1.9*1.9*2.8*1800*600*256)
|
|---|
| 187 | plot2d
|
|---|
| 188 | # plot2d lsspk x val*$VOL 1 'logx logy nsta xylimits=0.01,2.,10.,1e4 cpts marker=box,5 gold'
|
|---|
| 189 | plot2d lsspklobewn x val*$VOL 1 'same nsta cpts marker=box,5 siennared'
|
|---|
| 190 | plot2d subpkcorlobe x val*$VOL 1 'same nsta cpts marker=box,5 red'
|
|---|
| 191 | plot2d subpklobe x val*$VOL 1 'same nsta cpts marker=box,5 blueviolet'
|
|---|
| 192 | plot2d subpknolss x val*$VOL 1 'same nsta cpts marker=box,5 green'
|
|---|
| 193 |
|
|---|
| 194 | # settitle 'Recovered Pk[LSS] In=LSS+(GSM) (D=50 m)' ' ' 'font=helvetica,bold,18'
|
|---|
| 195 | settitle 'Recovered Pk[LSS] In=LSS+(Haslam+North20cm) (D=50 m)' ' ' 'font=helvetica,bold,18'
|
|---|
| 196 | setaxelabels 'k (Mpc^-1) h=0.7' 'P(k) (mK^2 Mpc^3)' 'font=helvetica,bolditalic,16'
|
|---|
| 197 | set lines ( 'Pk[LSS]' 'Pk[LSS*lobe+noise]' 'Pk[ExtractedLSS]' 'Pk[ExtLSS,NoBeamCor]' 'Pk[residual,NoLSS]' )
|
|---|
| 198 | set cols ( gold siennared red blueviolet green )
|
|---|
| 199 | textdrawer lines cols 'font=helvetica,bold,16 frame'
|
|---|
| 200 |
|
|---|
| 201 | plot2d fgndpk x val*$VOL 1 'logx logy xylimits=0.01,1.,1.,1e10 nsta cpts marker=box,5 black'
|
|---|
| 202 | plot2d fgndpkflobe x val*$VOL 1 ' nsta cpts marker=circle,5 navyblue same'
|
|---|
| 203 | plot2d fgndpklobe x val*$VOL 1 ' nsta cpts marker=circle,5 blue same'
|
|---|
| 204 |
|
|---|
| 205 | plot2d lsspk x val*$VOL 1 ' nsta cpts marker=box,5 red same'
|
|---|
| 206 | plot2d lsspkflobe x val*$VOL 1 ' nsta cpts marker=circle,5 orange same'
|
|---|
| 207 | plot2d lsspklobe x val*$VOL 1 ' nsta cpts marker=circle,5 yellow same'
|
|---|
| 208 |
|
|---|
| 209 | c++exec \
|
|---|
| 210 | Histo lsspkratioA = subpkcorlobe/lsspk; KeepObj(lsspkratioA); \
|
|---|
| 211 | Histo lsspkratioB = subpkcorlobe/lsspkflobe; KeepObj(lsspkratioB);
|
|---|
| 212 |
|
|---|
| 213 | plot2d lsspk x val*$VOL 1 'logx logy nsta xylimits=0.01,2.,10.,1e4 cpts marker=box,5 gold'
|
|---|
| 214 | plot2d lsspklobewn x val*$VOL 1 'same nsta cpts marker=box,5 red'
|
|---|
| 215 | plot2d subpknolss x val*$VOL 1 'same nsta cpts marker=box,5 green'
|
|---|
| 216 | plot2d subpknolssnocor x val*$VOL 1 'same nsta cpts marker=box,5 magenta'
|
|---|
| 217 | setaxelabels 'k (Mpc^-1) h=0.7' 'P(k) (mK^2 Mpc^3)' 'font=helvetica,bolditalic,16'
|
|---|
| 218 | settitle 'Recovered Pk[LSS] and residual systematics' ' ' 'font=helvetica,bold,18'
|
|---|
| 219 | set lines ( 'Pk[LSS]' 'Pk[LSS*lobe+noise]' 'Pk[residual,NoLSS]' 'Pk[residual,NoLSS,NoBeamCorrection]' )
|
|---|
| 220 | set cols ( gold red green magenta )
|
|---|
| 221 | textdrawer lines cols 'font=helvetica,bold,16 frame'
|
|---|