CosmoCoffee Forum Index CosmoCoffee

 
 FAQFAQ   SearchSearch  MemberlistSmartFeed   MemberlistMemberlist    RegisterRegister 
   ProfileProfile   Log inLog in 
Arxiv New Filter | Bookmarks & clubs | Arxiv ref/author:

COSMOMC: all_l_exact data format
Goto page 1, 2  Next  
Post new topic   Reply to topic    CosmoCoffee Forum Index -> Computers and software
View previous topic :: View next topic  
Author Message
Tor Helge Huse



Joined: 15 Dec 2004
Posts: 6
Affiliation: Astrofysisk Institutt, Oslo

PostPosted: April 18 2005  Reply with quote

I'm trying to use the all_l_exact data format in cosmomc with some simulated data, but I keep running into computational errors when doing so. As I can't find any documentation on what the all_l_exact data format should be, I'll ask here: What does the all_l_exact data format actually look like?

regards
Tor Helge Huse

PS: I want to first use CAMB to generate some data and then introduce some noise through my own routine, so CAMB data is my starting point.
Back to top
View user's profile  
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: April 18 2005  Reply with quote

The all_l_exact format uses the exact full-sky likelihood function given in e.g. Eq. 5 of astro-ph/0502469.

A .dataset file for an all_l_exact simulation would look something like

Code:
name = Planck
has_pol = T
all_l_exact=T
all_l_file = data/Planck.dat
all_l_lmax=2000


Then you just need to make the Planck.dat file with the mock data. The columns in this file are

l CTT (CTE CEE [CBB]) NT (NP) fskyeff

The brackets are included if has_pol = T, and CBB if you have compiled cosmomc with num_cls=4.

CTT (and the polarized ones similarly) are raw Cl + noise (no l(l+1) or 2π factors). i.e. they increase on scales below the beam size due to the larger noise there.

Here's some simple code to make a Planck.dat file from a Cl file output by CAMB (in μK2):

Code:
!Simulate CMB C_l data
!Simple version assuming isotropic noise and full sky
!no lensing, no tensors (check input Cl file has only 4 columns!)
!all parameters hard wired in code
!Data points are max likelihood values
!AL: Sept 2004

program CMB_cl_sim
 use AMLutils
 implicit none

!parameters
 character(LEN=*), parameter :: cl_infile = 'MAPcls.dat'
 character(LEN=*), parameter :: out_name = 'planck.dat'
 logical :: DoPol = .true., DoTemp = .true.
 real :: NoiseVar = 3d-4   !Pessimistic Planck
       !in MicroK^2
 real, parameter :: fwhm_arcmin = 7. ! 13.  !WMAP ~ 13 (0.22 deg min)
 real, parameter :: ENoiseFac = 2 !factor more noise on the polarization
 integer, parameter :: lmax = 2250


!Other stuff
 integer, parameter ::  C_T = 1, C_E = 2, C_B = 3, C_C = 4
 
 real :: fwhm_deg = fwhm_arcmin/60
 real :: Cl(lmax,4), NoiseCl(lmax)
 real amp, xlc, sigma2
 integer l,ll
 real T, E, TE,B, scal


print *,'lmax = ',lmax
print *,'fwhm_arcmin = ', fwhm_arcmin
print *,'Noise = ',NoiseVar

   xlc= 180*sqrt(8.*log(2.))/3.14159
   sigma2 = (fwhm_deg/xlc)**2
   do l=2, lmax
     NoiseCl(l) = NoiseVar*exp(l*(l+1)*sigma2)
   end do


Cl=0
B=0
open(unit=1,file=Cl_infile,form='formatted',status='old')
do ll=2,lmax
  read (1,*) l, T, E, TE
  if (l/=ll) stop 'cl must start at l=2'
  scal=twopi/(l*(l+1))
  Cl(l,1) = T*scal
  Cl(l,2) = E*scal
  Cl(l,3) = B*scal
  Cl(l,4) = TE*scal
end do
close(1)

call CreateTxtFile(out_name,34)

 do l=2, lmax

      write (34,'(1I5,6E15.5)') l,(Cl(l,1)+NoiseCl(l)) , &
                                (Cl(l,C_C)), &
                                (Cl(l,C_E)+ENoiseFac*NoiseCl(l)), &
                            !      (Cl(l,C_B)+ENoiseFac*Noise%Cl(l,1))*1e6, &
                                    NoiseCl(l)*1e6,ENoiseFac*NoiseCl(l), 1.
 end do
 close(34)
end program CMB_cl_sim


If you want to use a proper simulation rather than taking the data to be the exact model, you should use 1/(2l+1) m |alm|2 in the output, where alm is a Gaussian generated with variance Cl + Nl, or general alm from the theory Cl only and then add your own noise.


Last edited by Antony Lewis on April 19 2005; edited 1 time in total
Back to top
View user's profile [ Hidden ] Visit poster's website
Tor Helge Huse



Joined: 15 Dec 2004
Posts: 6
Affiliation: Astrofysisk Institutt, Oslo

PostPosted: April 18 2005  Reply with quote

Thank you for the information.
I'm trying to include tensors and when I do I get a very high Cl_BB term (~300 for l=20) for the first 100 or so l values, then dropping. Is this as it should be? The NT term is also very large (same order of magnitude as Cl_BB).
When running with the generated data the likelihood values get very large (xE+9) and it does not seem to get anywhere.

Just to clarify: should COBE_normalize be set to true in CAMB?

regards

Tor Helge
Back to top
View user's profile  
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: April 18 2005  Reply with quote

I have not tried it with tensors. Obviosuly the relative size of the noise depends on your experiment parameter, though for Planck it is correct that for sensible tensor amplitudes the signal will be noise dominated.

I wouldn't recommend using COBE_normalize, though you can as long as you make sure the final Cl is in μK2.

The noise free BB Cl should be about 0.001 the size of the temperature ones, so maybe you have a bug (not using the right column of the totCls output file or something: BB is column 4 of an output file with BB in).
Back to top
View user's profile [ Hidden ] Visit poster's website
Tor Helge Huse



Joined: 15 Dec 2004
Posts: 6
Affiliation: Astrofysisk Institutt, Oslo

PostPosted: April 19 2005  Reply with quote

Quote:

Code:

 do l=2, lmax

      write (34,'(1I5,6E15.5)') l,(Cl(l,1)+NoiseCl(l)) , &
                                (Cl(l,C_C)), &
                                (Cl(l,C_E)+ENoiseFac*NoiseCl(l)), &
                            !      (Cl(l,C_B)+ENoiseFac*Noise%Cl(l,1))*1e6, &
                                    NoiseCl(l)*1e6,ENoiseFac*NoiseCl(l), 1.
 end do
 close(34)
end program CMB_cl_sim



I'm a bit confused when it comes to the above. Why is the Cl(B) term and NoiseCl(l) term multiplied by 1e6? This blows the numbers up a lot, and especially the noise term, and as far as I can that term is used unchanged in cosmomc.

And the line
Code:

 !      (Cl(l,C_B)+ENoiseFac*Noise%Cl(l,1))*1e6, &

is for including tensors right? Is the term
Code:
Noise%Cl(l,1)
correct? I recieve a compilation error as it is not defined anywhere. Should it be
Code:
NoiseCl(l)
instead?

Tor Helge

Edit:
Quote:

I wouldn't recommend using COBE_normalize, though you can as long as you make sure the final Cl is in μK2.

Now running with a fresh install of CAMB (april05), nothing changed except lmax, kmax and ns (0.96), nt(−0.04) and amp_ratio(0.248).
Quote:

The noise free BB Cl should be about 0.001 the size of the temperature ones, so maybe you have a bug (not using the right column of the totCls output file or something: BB is column 4 of an output file with BB in).

The Cl(BB) is approx. 0.00001*Cl(TT) in the output from CAMB (column 4 vs column 2, l=48).
After running the data through your routine Cl(BB) is approx. 0.0001*Cl(TT)(l=48).
Back to top
View user's profile  
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: April 19 2005  Reply with quote

Yes, sorry, I missed deleting those 1e6 factors when converting from my code (which originally used Healpix-normalized output in mK2). Simiarly with the incompatible name. Obviously you need to replace CreateTxtFile as well unless you are using my utils library.

Should be

Code:
!Simulate CMB C_l data
!Simple version assuming isotropic noise and full sky
!no lensing, no tensors (check input Cl file has only 4 columns!)
!all parameters hard wired in code
!Data points are max likelihood values
!AL: Sept 2004

program CMB_cl_sim
 use AMLutils
 implicit none

!parameters
 character(LEN=*), parameter :: cl_infile = 'MAPcls.dat'
 character(LEN=*), parameter :: out_name = 'planck.dat'
 logical :: DoPol = .true., DoTemp = .true.
 real :: NoiseVar = 3d-4   !Pessimistic Planck
       !in MicroK^2
 real, parameter :: fwhm_arcmin = 7. ! 13.  !WMAP ~ 13 (0.22 deg min)
 real, parameter :: ENoiseFac = 2 !factor more noise on the polarization
 integer, parameter :: lmax = 2250


!Other stuff
 integer, parameter ::  C_T = 1, C_E = 2, C_B = 3, C_C = 4
 
 real :: fwhm_deg = fwhm_arcmin/60
 real :: Cl(lmax,4), NoiseCl(lmax)
 real amp, xlc, sigma2
 integer l,ll
 real T, E, TE,B, scal


print *,'lmax = ',lmax
print *,'fwhm_arcmin = ', fwhm_arcmin
print *,'Noise = ',NoiseVar

   xlc= 180*sqrt(8.*log(2.))/3.14159
   sigma2 = (fwhm_deg/xlc)**2
   do l=2, lmax
     NoiseCl(l) = NoiseVar*exp(l*(l+1)*sigma2)
   end do


Cl=0
B=0
open(unit=1,file=Cl_infile,form='formatted',status='old')
do ll=2,lmax
  read (1,*) l, T, E, TE
  if (l/=ll) stop 'cl must start at l=2'
  scal=twopi/(l*(l+1))
  Cl(l,1) = T*scal
  Cl(l,2) = E*scal
  Cl(l,3) = B*scal
  Cl(l,4) = TE*scal
end do
close(1)

call CreateTxtFile(out_name,34)

 do l=2, lmax

      write (34,'(1I5,6E15.5)') l,(Cl(l,1)+NoiseCl(l)) , &
                                (Cl(l,C_C)), &
                                (Cl(l,C_E)+ENoiseFac*NoiseCl(l)), &
                            !      (Cl(l,C_B)+ENoiseFac*NoiseCl(l)), &
                                    NoiseCl(l),ENoiseFac*NoiseCl(l), 1.
 end do
 close(34)
end program CMB_cl_sim
Back to top
View user's profile [ Hidden ] Visit poster's website
Tor Helge Huse



Joined: 15 Dec 2004
Posts: 6
Affiliation: Astrofysisk Institutt, Oslo

PostPosted: April 19 2005  Reply with quote

Ok, thank you. That makes a whole lot more sense. Testing it now (fingers crossed).

Tor Helge
Back to top
View user's profile  
Samuel Leach



Joined: 15 Oct 2004
Posts: 16
Affiliation: SISSA, Trieste

PostPosted: May 11 2005  Reply with quote

Dear Antony,

I've been following your discussion here and have been trying out your all_l_exact data format with CosmoMC which of course works rather superbly. Just a couple of related questions here that I'd like you to help me think through:

- I can see that using the exact C_l as the data is a Gaussian realisation, and a rather special one which can be useful for testing the mechanics of the data analysis pipeline for intrinsic bias in parameter reconstruction and in the determination of derived parameters, for instance. Do you know whether this will always lead to the right assessment for the posterior volume ? Of course I should just implement the Gaussian realisation version of the code to experiment a bit.

- What about your Planck noise model ? Can we set up a small database on CosmoCoffee which contains the current broad expectation for Planck and other future experiments, in terms of noise variance and beam width?

- What would be the standard way to cite all the useful information, codes etc from CosmoCofeee ? Did you discuss this together already? These aren't exactly Private Communications.

Thanks always for your time and energies Antony.
Back to top
View user's profile [ Hidden ] Visit poster's website
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: May 11 2005  Reply with quote

Samuel Leach wrote:
- I can see that using the exact Cl as the data is a Gaussian realisation, and a rather special one which can be useful for testing the mechanics of the data analysis pipeline for intrinsic bias in parameter reconstruction and in the determination of derived parameters, for instance. Do you know whether this will always lead to the right assessment for the posterior volume ? Of course I should just implement the Gaussian realisation version of the code to experiment a bit.


This is an interesting question, and I don't really know the answer. It probably depends on what you think the "right" answer is for forecasting work when you don't know what the exact realisation will be.

My experience for the CMB is that the error bars come out very similar to those you get from a typical realisation (with the advantage that the peak should be exactly at the input parameter values). However, I've not been able to prove any relation to any exact Bayesian forcasting protocol – I'd be interested to know any thoughts people have.

Note that that actual quadrupole appears to be low, and hence not very typical. Error bars can depend slightly on the value of the quadrupole realisation, which is an example where using best-fit theory value may not be realistic for future experiment error bars (given what we already know about the realisation from WMAP).

Samuel Leach wrote:

- What about your Planck noise model ? Can we set up a small database on CosmoCoffee which contains the current broad expectation for Planck and other future experiments, in terms of noise variance and beam width?


Sarah started work on a Wiki giving details of future experiments etc. This would be useful, but of course someone has to put in the work to input all the data. (though the idea with a Wiki of course is that anyone can do it)

Samuel Leach wrote:

- What would be the standard way to cite all the useful information, codes etc from CosmoCofeee ? Did you discuss this together already? These aren't exactly Private Communications.


I've referred to plots on CosmoCoffee - just cited the URL. Otherwise I guess mention CosmoCoffee in the Acknowledgements - but we haven't really discussed this.
Back to top
View user's profile [ Hidden ] Visit poster's website
Pascal Vaudrevange



Joined: 26 Mar 2006
Posts: 50
Affiliation: DESY

PostPosted: May 23 2007  Reply with quote

Reviving this old thread.. I would like to create some mock datasets with WMAP3 accuracy. What would the values for noise_var be? Thanks
Back to top
View user's profile [ Hidden ] Visit poster's website
Pascal Vaudrevange



Joined: 26 Mar 2006
Posts: 50
Affiliation: DESY

PostPosted: May 24 2007  Reply with quote

I guess I can answer the question myself. I found this paper by Antony

http://arxiv.org/abs/astro-ph/0502469

which uses NlTT=0.03μ K2 . (Of course, this is only a crude approximation as the beam is not Gaussian..)
Back to top
View user's profile [ Hidden ] Visit poster's website
Jochen Weller



Joined: 24 Sep 2004
Posts: 43
Affiliation: Ludwig-Maximilians-University Munich

PostPosted: March 08 2008  Reply with quote

I have a question regarding the input format of the exact Cl mock data.
In cmbdata.F90 you write:

Code:
  !File format:
     !l C_TT (C_TE C_EE [C_BB]) N_T (N_P) fsky_eff
     !Must have num_cls set to correct value for file
         do l = 2, aset%all_l_lmax
            read (file_unit, *, end=100, err=100) idum,inobs(1:ncls), aset%all_l_noise(l,:),
            aset%all_l_fsky(l)
           if (idum /= l) stop 'Error reading all_l_file'
            !set BB to pure noise if not in file
          if (aset%has_pol .and. ncls < num_cls) inobs(num_cls) = aset%all_l_noise(l,2)
          aset%all_l_obs(l,:) = inobs(1:num_cls)
        end do



However in the example code to create the mock data you write:
Code:

     do l=2, lmax

      write (34,'(1I5,6E15.5)') l,(Cl(l,1)+NoiseCl(l)) , &
                                (Cl(l,C_C)), &
                                (Cl(l,C_E)+ENoiseFac*NoiseCl(l)), &
                            !      (Cl(l,C_B)+ENoiseFac*NoiseCl(l)), &
                                    NoiseCl(l),ENoiseFac*NoiseCl(l), 1.
       end do



Now I am a bit puzzled by why the second column has to be written Cl(l,1)+NoiseCl(l)
I would have expected to have Cl(l,C_Temp) properly normalized to be in there, from the above read statement?
Back to top
View user's profile [ Hidden ] Visit poster's website
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: March 08 2008  Reply with quote

I'm not sure I get the question, but all Cl here are the bare ones without factors of l(l+1)/2π.

C_Temp=1.
Back to top
View user's profile [ Hidden ] Visit poster's website
Jochen Weller



Joined: 24 Sep 2004
Posts: 43
Affiliation: Ludwig-Maximilians-University Munich

PostPosted: March 08 2008  Reply with quote

What I meant it looks like in the 2nd column you are writing Cl+noise, but you just read Cl
Back to top
View user's profile [ Hidden ] Visit poster's website
Antony Lewis



Joined: 23 Sep 2004
Posts: 1237
Affiliation: University of Sussex

PostPosted: March 08 2008  Reply with quote

Yes, the Cl in the data file are the 'observed' Cl and include noise. The exact likelihood is a function of only Cl + Nl for full sky isotropic noise.
Back to top
View user's profile [ Hidden ] Visit poster's website
Display posts from previous:   
Post new topic   Reply to topic    CosmoCoffee Forum Index -> Computers and software All times are GMT + 5 Hours
Goto page 1, 2  Next
Page 1 of 2

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group. Sponsored by WordWeb online dictionary and dictionary software.