CosmicFish  Reference documentation for version 1.0
Looking into future Cosmology
005_init_from_file.f90
Go to the documentation of this file.
1 !----------------------------------------------------------------------------------------
2 !
3 ! This file is part of CosmicFish.
4 !
5 ! Copyright (C) 2015-2016 by the CosmicFish authors
6 !
7 ! The CosmicFish code is free software;
8 ! You can use it, redistribute it, and/or modify it under the terms
9 ! of the GNU General Public License as published by the Free Software Foundation;
10 ! either version 3 of the License, or (at your option) any later version.
11 ! The full text of the license can be found in the file LICENSE at
12 ! the top level of the CosmicFish distribution.
13 !
14 !----------------------------------------------------------------------------------------
15 
19 
20 !----------------------------------------------------------------------------------------
23 
25 
26 module init_from_file
27 
28  use precision
29  use cosmicfish_types
30  use inifile
31  use modelparams
32 
33  implicit none
34 
35  private
36 
37  public init_cosmicfish_from_file
38 
39 contains
40 
41  ! ---------------------------------------------------------------------------------------------
43  subroutine init_cosmicfish_from_file( P, FP, filename, param_out_name )
44 
45  use constants, only : f_21cm, cobe_cmbtemp
46  use camb
47  use cambmain, only : alens
48  use bispectrum
49 
50  implicit none
51 
52  Type(cambparams) :: P
53  Type(cosmicfish_params) :: FP
54  character(len=*), intent(in) :: filename
55  character(len=*), intent(in), optional :: param_out_name
56 
57  logical :: bad
58  logical :: DoCounts = .false.
59  integer :: i, ind
60  type(tredwin), pointer :: redwin
61  character(LEN=Ini_max_string_len) numstr, S, outroot, version_check
62  real(dl) :: nmassive
63  character(LEN=:), allocatable :: ParamDir
64 
65  !MMmod
66  character(LEN=Ini_max_string_len) red_ind
67 
68  ! open the parameter file:
69  call ini_open(filename, 1, bad, .false.)
70  if (bad) stop 'Error opening parameter file'
71  ini_fail_on_not_found = .false.
72  ! output root:
73  outroot = ini_read_string('output_root')
74  if (outroot /= '') outroot = trim(outroot) // '_'
75  ! set CAMB parameters to defaul at the beginning:
76  call camb_setdefparams(p)
77  ! high l template:
78  if (ini_haskey('highL_unlensed_cl_template')) then
79  highl_unlensed_cl_template= ini_read_string_default('highL_unlensed_cl_template',highl_unlensed_cl_template)
80  else
81  ! assume that the file highL_unlensed_cl_template is in the parameter folder:
82  i = scan( filename, '/' , .true. )
83  paramdir = ''
84  if ( i .gt. 1 ) then
85  paramdir = filename(1:i)
86  endif
87  highl_unlensed_cl_template = concat(paramdir, highl_unlensed_cl_template)
88  end if
89  ! read what to do:
90  p%WantScalars = ini_read_logical('get_scalar_cls',.true.)
91  p%WantVectors = ini_read_logical('get_vector_cls',.false.)
92  p%WantTensors = ini_read_logical('get_tensor_cls',.false.)
93 
94  p%Want_CMB = ini_read_logical('want_CMB',.true.)
95  p%Want_CMB_lensing = ini_read_logical('want_CMB_lensing',.true.)
96  p%Want_CMB_lensing = p%Want_CMB .or. p%Want_CMB_lensing
97 
98  p%WantCls= p%WantScalars .or. p%WantTensors .or. p%WantVectors
99  p%PK_WantTransfer=ini_read_logical('get_transfer')
100  ! read the window functions:
101  if (p%WantScalars) then
102  num_redshiftwindows = ini_read_int('num_redshiftwindows',0)
103  else
104  num_redshiftwindows = 0
105  end if
106  limber_windows = ini_read_logical('limber_windows',limber_windows)
107  if (limber_windows) limber_phiphi = ini_read_int('limber_phiphi',limber_phiphi)
108  if (num_redshiftwindows>0) then
109  doredshiftlensing = ini_read_logical('DoRedshiftLensing',.false.)
110  kmax_boost = ini_read_double('Kmax_Boost',kmax_boost)
111  end if
112  do21cm = ini_read_logical('Do21cm', .false.)
113  num_extra_redshiftwindows = 0
114 
115  ! initialize the window functions:
116  do i=1, num_redshiftwindows
117  redwin => redshift_w(i)
118  call initredshiftwindow(redwin)
119  write (numstr,*) i
120  numstr=adjustl(numstr)
121  redwin%Redshift = ini_read_double('redshift('//trim(numstr)//')')
122  s = ini_read_string('redshift_kind('//trim(numstr)//')')
123  if (s=='21cm') then
124  redwin%kind = window_21cm
125  elseif (s=='counts') then
126  redwin%kind = window_counts
127  elseif (s=='lensing') then
128  redwin%kind = window_lensing
129  else
130  write (*,*) i, 'Error: unknown type of window '//trim(s)
131  stop
132  end if
133  redwin%a = 1/(1+redwin%Redshift)
134  if (redwin%kind /= window_21cm) then
135  redwin%sigma = ini_read_double('redshift_sigma('//trim(numstr)//')')
136  redwin%sigma_z = redwin%sigma
137  else
138  do21cm = .true.
139  redwin%sigma = ini_read_double('redshift_sigma_Mhz('//trim(numstr)//')')
140  if (redwin%sigma < 0.003) then
141  write(*,*) 'WARNING:Window very narrow.'
142  write(*,*) ' --> use transfer functions and transfer_21cm_cl =T ?'
143  end if
144  !with 21cm widths are in Mhz, make dimensionless scale factor
145  redwin%sigma = redwin%sigma/(f_21cm/1e6)
146  redwin%sigma_z = redwin%sigma*(1+redwin%RedShift)**2
147  write(*,*) i,'delta_z = ', redwin%sigma_z
148  end if
149  if (redwin%kind == window_counts) then
150  docounts = .true.
151  redwin%bias = ini_read_double('redshift_bias('//trim(numstr)//')')
152  redwin%dlog10Ndm = ini_read_double('redshift_dlog10Ndm('//trim(numstr)//')',0.d0)
153  if (doredshiftlensing) then
154  num_extra_redshiftwindows=num_extra_redshiftwindows+1
155  redwin%mag_index = num_extra_redshiftwindows
156  end if
157  end if
158  end do
159 
160  if (do21cm) then
161  print*, 'Not yet implemented'
162  stop
163  line_basic = ini_read_logical('line_basic')
164  line_distortions = ini_read_logical('line_distortions')
165  line_extra = ini_read_logical('line_extra')
166 
167  line_phot_dipole = ini_read_logical('line_phot_dipole')
168  line_phot_quadrupole = ini_read_logical('line_phot_quadrupole')
169  line_reionization = ini_read_logical('line_reionization')
170 
171  use_mk = ini_read_logical('use_mK')
172  if (debugmsgs) then
173  write (*,*) 'Doing 21cm'
174  write (*,*) 'dipole = ',line_phot_dipole, ' quadrupole =', line_phot_quadrupole
175  end if
176  else
177  line_extra = .false.
178  end if
179 
180  if (docounts) then
181  counts_density = ini_read_logical('counts_density')
182  counts_redshift = ini_read_logical('counts_redshift')
183  counts_radial = ini_read_logical('counts_radial')
184  counts_evolve = ini_read_logical('counts_evolve')
185  counts_timedelay = ini_read_logical('counts_timedelay')
186  counts_isw = ini_read_logical('counts_ISW')
187  counts_potential = ini_read_logical('counts_potential')
188  counts_velocity = ini_read_logical('counts_velocity')
189  end if
190 
191  p%OutputNormalization=outnone
192 
193  accuracyboost = ini_read_double('accuracy_boost',accuracyboost)
194  laccuracyboost = ini_read_real('l_accuracy_boost',laccuracyboost)
195  highaccuracydefault = ini_read_logical('high_accuracy_default',highaccuracydefault)
196 
197  p%NonLinear = ini_read_int('do_nonlinear',nonlinear_none)
198 
199  evolve_delta_xe = ini_read_logical('evolve_delta_xe', .false.)
200 
201  p%DoLensing = .false.
202  if (p%WantCls) then
203  if (p%WantScalars .or. p%WantVectors) then
204  p%Max_l = ini_read_int('l_max_scalar')
205  p%Max_eta_k = ini_read_double('k_eta_max_scalar',p%Max_l*2._dl)
206  if (p%WantScalars) then
207  p%DoLensing = ini_read_logical('do_lensing',.false.)
208  if (p%DoLensing) lensing_method = ini_read_int('lensing_method',1)
209  end if
210  if (p%WantVectors) then
211  if (p%WantScalars .or. p%WantTensors) stop 'Must generate vector modes on their own'
212  i = ini_read_int('vector_mode')
213  if (i==0) then
214  vec_sig0 = 1
215  magnetic = 0
216  else if (i==1) then
217  magnetic = -1
218  vec_sig0 = 0
219  else
220  stop 'vector_mode must be 0 (regular) or 1 (magnetic)'
221  end if
222  end if
223  end if
224 
225  if (p%WantTensors) then
226  p%Max_l_tensor = ini_read_int('l_max_tensor')
227  p%Max_eta_k_tensor = ini_read_double('k_eta_max_tensor',max(500._dl,p%Max_l_tensor*2._dl))
228  end if
229  endif
230 
231  p%window_kmax_boost = ini_read_double('window_kmax_boost', 1._dl)
232 
233 #ifdef COSMICFISH_CAMB
234  ! Read cosmological parameters.
235  call darkenergy_readparams( p, defini )
236 #endif
237 
238  p%h0 = ini_read_double('hubble')
239 
240  if (ini_read_logical('use_physical',.false.)) then
241  p%omegab = ini_read_double('ombh2')/(p%H0/100)**2
242  p%omegac = ini_read_double('omch2')/(p%H0/100)**2
243  p%omegan = ini_read_double('omnuh2')/(p%H0/100)**2
244  p%omegav = 1- ini_read_double('omk') - p%omegab-p%omegac - p%omegan
245  else
246  p%omegab = ini_read_double('omega_baryon')
247  p%omegac = ini_read_double('omega_cdm')
248  p%omegav = ini_read_double('omega_lambda')
249  p%omegan = ini_read_double('omega_neutrino')
250  end if
251 
252  p%tcmb = ini_read_double('temp_cmb',cobe_cmbtemp)
253 
254 #ifdef COSMICFISH_EFTCAMB
255  ! Read EFTCAMB parameters.
256 
257  ! 1) Initialization of EFTCAMB flags.
258 
259  p%EFTflag = ini_read_int('EFTflag',0)
260 
261  p%PureEFTmodelOmega = ini_read_int('PureEFTmodelOmega',0)
262  p%PureEFTmodelGamma1 = ini_read_int('PureEFTmodelGamma1',0)
263  p%PureEFTmodelGamma2 = ini_read_int('PureEFTmodelGamma2',0)
264  p%PureEFTmodelGamma3 = ini_read_int('PureEFTmodelGamma3',0)
265  p%PureEFTmodelGamma4 = ini_read_int('PureEFTmodelGamma4',0)
266  p%PureEFTmodelGamma5 = ini_read_int('PureEFTmodelGamma5',0)
267  p%PureEFTmodelGamma6 = ini_read_int('PureEFTmodelGamma6',0)
268 
269  p%DesignerEFTmodel = ini_read_int('DesignerEFTmodel',1)
270  p%AltParEFTmodel = ini_read_int('AltParEFTmodel',1)
271  p%FullMappingEFTmodel = ini_read_int('FullMappingEFTmodel',1)
272 
273  ! 2) Initialization of EFTCAMB model properties flags.
274 
275  ! read the DE eos model selection flag:
276  p%EFTwDE = ini_read_int('EFTwDE',0)
277  ! read pure EFT Horndeski model selection flag:
278  p%PureEFTHorndeski = ini_read_logical('PureEFTHorndeski',.false.)
279  ! read RPH model selection flags:
280  p%RPHmassPmodel = ini_read_int('RPHmassPmodel',0)
281  p%RPHkineticitymodel = ini_read_int('RPHkineticitymodel',0)
282  p%RPHbraidingmodel = ini_read_int('RPHbraidingmodel',0)
283  p%RPHtensormodel = ini_read_int('RPHtensormodel',0)
284  ! read the Horava Solar System Free flag:
285  p%HoravaSolarSystem = ini_read_logical('HoravaSolarSystem',.false.)
286 
287  ! 3) Initialization of EFTCAMB stability flags:
288 
289  p%EFT_mathematical_stability = ini_read_logical('EFT_mathematical_stability',.true.)
290  p%EFT_physical_stability = ini_read_logical('EFT_physical_stability',.true.)
291  p%EFTAdditionalPriors = ini_read_logical('EFTAdditionalPriors',.true.)
292  p%MinkowskyPriors = ini_read_logical('MinkowskyPriors',.true.)
293 
294  ! 4) Initialization of EFTCAMB model parameters.
295 
296  ! read the DE eos parameters:
297  p%EFTw0 = ini_read_double('EFTw0',-1._dl)
298  p%EFTwa = ini_read_double('EFTwa',0._dl)
299  p%EFTwn = ini_read_double('EFTwn',2._dl)
300  p%EFTwat = ini_read_double('EFTwat',1._dl)
301  p%EFtw2 = ini_read_double('EFtw2',0._dl)
302  p%EFTw3 = ini_read_double('EFTw3',0._dl)
303  ! read pure EFT parameters:
304  p%EFTOmega0 = ini_read_double('EFTOmega0', 0.0_dl)
305  p%EFTOmegaExp = ini_read_double('EFTOmegaExp', 0.0_dl)
306  p%EFTGamma10 = ini_read_double('EFTGamma10', 0.0_dl)
307  p%EFTGamma1Exp = ini_read_double('EFTGamma1Exp', 0.0_dl)
308  p%EFTGamma20 = ini_read_double('EFTGamma20', 0.0_dl)
309  p%EFTGamma2Exp = ini_read_double('EFTGamma2Exp', 0.0_dl)
310  p%EFTGamma30 = ini_read_double('EFTGamma30', 0.0_dl)
311  p%EFTGamma3Exp = ini_read_double('EFTGamma3Exp', 0.0_dl)
312  p%EFTGamma40 = ini_read_double('EFTGamma40', 0.0_dl)
313  p%EFTGamma4Exp = ini_read_double('EFTGamma4Exp', 0.0_dl)
314  p%EFTGamma50 = ini_read_double('EFTGamma50', 0.0_dl)
315  p%EFTGamma5Exp = ini_read_double('EFTGamma5Exp', 0.0_dl)
316  p%EFTGamma60 = ini_read_double('EFTGamma60', 0.0_dl)
317  p%EFTGamma6Exp = ini_read_double('EFTGamma6Exp', 0.0_dl)
318  ! read f(R) parameters:
319  p%EFTB0 = ini_read_double('EFTB0', 0.0_dl)
320  ! read RPH parameters:
321  p%RPHmassP0 = ini_read_double('RPHmassP0', 0.0_dl)
322  p%RPHmassPexp = ini_read_double('RPHmassPexp', 0.0_dl)
323  p%RPHkineticity0 = ini_read_double('RPHkineticity0', 0.0_dl)
324  p%RPHkineticityexp = ini_read_double('RPHkineticityexp', 0.0_dl)
325  p%RPHbraiding0 = ini_read_double('RPHbraiding0', 0.0_dl)
326  p%RPHbraidingexp = ini_read_double('RPHbraidingexp', 0.0_dl)
327  p%RPHtensor0 = ini_read_double('RPHtensor0', 0.0_dl)
328  p%RPHtensorexp = ini_read_double('RPHtensorexp', 0.0_dl)
329  ! read Horava parameters:
330  p%Horava_xi = ini_read_double('Horava_xi', 0.0_dl)
331  p%Horava_lambda = ini_read_double('Horava_lambda', 0.0_dl)
332  p%Horava_eta = ini_read_double('Horava_eta', 0.0_dl)
333 #endif
334 
335 #ifdef COSMICFISH_MGCAMB
336  p%MGC_model = ini_read_int('model',0)
337  p%GRtrans = ini_read_double('GRtrans',0.d0)
338 
339  if ( p%MGC_model == 1 ) then
340  p%B1 = ini_read_double('B1',0.d0)
341  p%B2 = ini_read_double('B2',0.d0)
342  p%lambda1_2 = ini_read_double('lambda1_2',0.d0)
343  p%lambda2_2 = ini_read_double('lambda2_2',0.d0)
344  p%ss = ini_read_double('ss',0.d0)
345  else if ( p%MGC_model == 2 ) then
346  p%MGQfix = ini_read_double('MGQfix',1.d0)
347  p%MGRfix = ini_read_double('MGRfix',1.d0)
348  else if ( p%MGC_model == 3 ) then
349  p%Qnot = ini_read_double('Qnot',1.d0)
350  p%Rnot = ini_read_double('Rnot',1.d0)
351  p%sss = ini_read_double('sss',0.d0)
352  else if ( p%MGC_model == 4 ) then
353  p%B0 = ini_read_double('B0',0.d0)
354  p%B1 = 4.d0/3.d0
355  p%lambda1_2 = (p%B0*(299792458.d-3)**2)/(2.d0*p%H0**2)
356  p%B2 = 0.5d0
357  p%lambda2_2 = p%B1*p%lambda1_2
358  p%ss = 4.d0
359  else if ( p%MGC_model ==5 ) then
360  p%B0 = ini_read_double('B0',0.d0)
361  p%B1 = ini_read_double('beta1',0.d0)
362  p%lambda1_2 = (p%B0*(299792458.d-3)**2)/(2.d0*p%H0**2)
363  p%B2 = 2.d0/p%B1 -1.d0
364  p%lambda2_2 = p%B1*p%lambda1_2
365  p%ss = ini_read_double('s',0.d0)
366  else if ( p%MGC_model ==6 ) then
367  p%Linder_gamma = ini_read_double('Linder_gamma',0.d0)
368  else if ( p%MGC_model == 7 ) then
369  p%beta_star = ini_read_double('beta_star', 0.d0)
370  p%xi_star = ini_read_double('xi_star', 0.d0)
371  p%a_star = ini_read_double('a_star', 0.d0)
372  p%GRtrans = p%a_star
373  else if ( p%MGC_model == 8 ) then
374  p%beta0 = ini_read_double('beta0', 0.d0)
375  p%xi0 = ini_read_double('xi0', 0.d0)
376  p%DilR = ini_read_double('DilR', 0.d0)
377  p%DilS = ini_read_double('DilS', 0.d0)
378  else if ( p%MGC_model == 9 ) then
379  p%F_R0 = ini_read_double('F_R0', 0.d0)
380  p%FRn = ini_read_double('FRn', 0.d0)
381  p%beta0 = 1.d0/sqrt(6.d0)
382  else if ( p%MGC_model ==10 ) then
383  p%beta0 = ini_read_double('beta0', 0.d0)
384  p%A_2 = ini_read_double('A2',0.d0)
385  else if ( p%MGC_model /= 0 ) then
386  print*, '***MGCAMB: please choose a model***'
387  stop
388  end if
389 #endif
390 
391  p%yhe = ini_read_double('helium_fraction',0.24_dl)
392  p%Num_Nu_massless = ini_read_double('massless_neutrinos')
393 
394  ! massive neutrinos:
395  p%Nu_mass_eigenstates = ini_read_int('nu_mass_eigenstates',1)
396  if (p%Nu_mass_eigenstates > max_nu) stop 'too many mass eigenstates'
397 
398  numstr = ini_read_string('massive_neutrinos')
399  read(numstr, *) nmassive
400  if (abs(nmassive-nint(nmassive))>1e-6) stop 'massive_neutrinos should now be integer (or integer array)'
401  read(numstr,*) p%Nu_Mass_numbers(1:p%Nu_mass_eigenstates)
402 
403  p%Num_Nu_massive = sum(p%Nu_Mass_numbers(1:p%Nu_mass_eigenstates))
404 
405  if (p%Num_Nu_massive>0) then
406  p%share_delta_neff = ini_read_logical('share_delta_neff', .true.)
407  numstr = ini_read_string('nu_mass_degeneracies')
408  if (p%share_delta_neff) then
409  if (numstr/='') write (*,*) 'WARNING: nu_mass_degeneracies ignored when share_delta_neff'
410  else
411  if (numstr=='') stop 'must give degeneracies for each eigenstate if share_delta_neff=F'
412  read(numstr,*) p%Nu_mass_degeneracies(1:p%Nu_mass_eigenstates)
413  end if
414  numstr = ini_read_string('nu_mass_fractions')
415  if (numstr=='') then
416  if (p%Nu_mass_eigenstates >1) stop 'must give nu_mass_fractions for the eigenstates'
417  p%Nu_mass_fractions(1)=1
418  else
419  read(numstr,*) p%Nu_mass_fractions(1:p%Nu_mass_eigenstates)
420  end if
421  end if
422 
423  !JD 08/13 begin changes for nonlinear lensing of CMB + LSS compatibility
424  !P%Transfer%redshifts -> P%Transfer%PK_redshifts and P%Transfer%num_redshifts -> P%Transfer%PK_num_redshifts
425  !in the P%WantTransfer loop.
426  if (((p%NonLinear==nonlinear_lens .or. p%NonLinear==nonlinear_both) .and. (p%DoLensing .or. num_redshiftwindows>0)) &
427  .or. p%PK_WantTransfer) then
428  p%Transfer%high_precision= ini_read_logical('transfer_high_precision',.false.)
429  else
430  p%transfer%high_precision = .false.
431  endif
432 
433  if (p%PK_WantTransfer) then
434  p%WantTransfer = .true.
435  p%transfer%kmax = ini_read_double('transfer_kmax')
436  p%transfer%k_per_logint = ini_read_int('transfer_k_per_logint')
437  p%transfer%PK_num_redshifts = ini_read_int('transfer_num_redshifts')
438 
439  if (do21cm) transfer_21cm_cl = ini_read_logical('transfer_21cm_cl',.false.)
440  if (transfer_21cm_cl .and. p%transfer%kmax > 800) then
441  !Actually line widths are important at significantly larger scales too
442  write (*,*) 'WARNING: kmax very large. '
443  write(*,*) ' -- Neglected line width effects will dominate'
444  end if
445 
446  transfer_interp_matterpower = ini_read_logical('transfer_interp_matterpower ', transfer_interp_matterpower)
447  transfer_power_var = ini_read_int('transfer_power_var',transfer_power_var)
448  if (p%transfer%PK_num_redshifts > max_transfer_redshifts) stop 'Too many redshifts'
449  do i=1, p%transfer%PK_num_redshifts
450  p%transfer%PK_redshifts(i) = ini_read_double_array('transfer_redshift',i,0._dl)
451  end do
452  else
453  p%Transfer%PK_num_redshifts = 1
454  p%Transfer%PK_redshifts = 0
455  end if
456 
457  if ((p%NonLinear==nonlinear_lens .or. p%NonLinear==nonlinear_both) .and. &
458  (p%DoLensing .or. num_redshiftwindows > 0)) then
459  p%WantTransfer = .true.
460  call transfer_setfornonlinearlensing(p%Transfer)
461  end if
462 
463  call transfer_sortandindexredshifts(p%Transfer)
464  !JD 08/13 end changes
465 
466  p%transfer%kmax=p%transfer%kmax*(p%h0/100._dl)
467 
468  ini_fail_on_not_found = .false.
469 
470  debugparam = ini_read_double('DebugParam',debugparam)
471  alens = ini_read_double('Alens',alens)
472 
473  call reionization_readparams(p%Reion, defini)
474  call initialpower_readparams(p%InitPower, defini, p%WantTensors)
475  call recombination_readparams(p%Recomb, defini)
476  if (ini_haskey('recombination')) then
477  i = ini_read_int('recombination',1)
478  if (i/=1) stop 'recombination option deprecated'
479  end if
480 
481  call bispectrum_readparams(bispectrumparams, defini, outroot)
482 
483  if (p%WantScalars .or. p%WantTransfer) then
484  p%Scalar_initial_condition = ini_read_int('initial_condition',initial_adiabatic)
485  if (p%Scalar_initial_condition == initial_vector) then
486  p%InitialConditionVector=0
487  numstr = ini_read_string('initial_vector',.true.)
488  read (numstr,*) p%InitialConditionVector(1:initial_iso_neutrino_vel)
489  end if
490  if (p%Scalar_initial_condition/= initial_adiabatic) use_spline_template = .false.
491  end if
492 
493  if (p%WantScalars) then
494  has_cl_2d_array = .true.
495  end if
496 
497  ini_fail_on_not_found = .false.
498 
499  !optional parameters controlling the computation
500  p%AccuratePolarization = ini_read_logical('accurate_polarization',.true.)
501  p%AccurateReionization = ini_read_logical('accurate_reionization',.false.)
502  p%AccurateBB = ini_read_logical('accurate_BB',.false.)
503  p%DerivedParameters = ini_read_logical('derived_parameters',.true.)
504 
505  version_check = ini_read_string('version_check')
506  if (version_check == '') then
507  !tag the output used parameters .ini file with the version of CAMB being used now
508  call tnamevaluelist_add(defini%ReadValues, 'version_check', version)
509  else if (version_check /= version) then
510  write(*,*) 'WARNING: version_check does not match this CAMB version'
511  end if
512  !Mess here to fix typo with backwards compatibility
513  if (ini_haskey('do_late_rad_trunction')) then
514  dolateradtruncation = ini_read_logical('do_late_rad_trunction',.true.)
515  if (ini_haskey('do_late_rad_truncation')) stop 'check do_late_rad_xxxx'
516  else
517  dolateradtruncation = ini_read_logical('do_late_rad_truncation',.true.)
518  end if
519  dotensorneutrinos = ini_read_logical('do_tensor_neutrinos',dotensorneutrinos )
520  feedbacklevel = ini_read_int('feedback_level',feedbacklevel)
521 
522  p%MassiveNuMethod = ini_read_int('massive_nu_approx',nu_best)
523 
524  threadnum = ini_read_int('number_of_threads',threadnum)
525  use_spline_template = ini_read_logical('use_spline_template',use_spline_template)
526 
527  dotensorneutrinos = dotensorneutrinos .or. highaccuracydefault
528  if (do_bispectrum) then
529  lsampleboost = 50
530  else
531  lsampleboost = ini_read_double('l_sample_boost',lsampleboost)
532  end if
533 
534  ! From now on read cosmicfish parameters:
535 
536  fp%outroot = trim(outroot)
537 
538  ! initialize general stuff:
539  fp%adaptivity = ini_read_logical('adaptivity',.false.)
540  fp%cosmicfish_feedback = ini_read_int('cosmicfish_feedback',0)
541 
542  fp%cosmicfish_want_cls = ini_read_logical('cosmicfish_want_cls',.false.)
543  fp%cosmicfish_want_SN = ini_read_logical('cosmicfish_want_SN' ,.false.)
544  fp%cosmicfish_want_Mpk = ini_read_logical('cosmicfish_want_Mpk',.false.)
545  fp%cosmicfish_want_RD = ini_read_logical('cosmicfish_want_RD' ,.false.)
546  fp%cosmicfish_want_derived = ini_read_logical('cosmicfish_want_derived',.false.)
547 
548  fp%output_factor = ini_read_double('CMB_outputscale',1.d0)
549 
550  ! read parameter selection:
551  fp%fisher_par%want_ombh2 = ini_read_logical( 'param[ombh2]' ,.false. )
552  fp%fisher_par%want_omch2 = ini_read_logical( 'param[omch2]' ,.false. )
553  fp%fisher_par%want_omnuh2 = ini_read_logical( 'param[omnuh2]' ,.false. )
554  fp%fisher_par%want_hubble = ini_read_logical( 'param[hubble]' ,.false. )
555  fp%fisher_par%want_helium_fraction = ini_read_logical( 'param[Ye]' ,.false. )
556  fp%fisher_par%want_massless = ini_read_logical( 'param[nom_nu]' ,.false. )
557  fp%fisher_par%want_scalar_amp = ini_read_logical( 'param[As]' ,.false. )
558  fp%fisher_par%want_scalar_spectral_index = ini_read_logical( 'param[ns]' ,.false. )
559  fp%fisher_par%want_scalar_nrun = ini_read_logical( 'param[nsrun]' ,.false. )
560  fp%fisher_par%want_tensor_spectral_index = ini_read_logical( 'param[nt]' ,.false. )
561  fp%fisher_par%want_initial_ratio = ini_read_logical( 'param[r]' ,.false. )
562  fp%fisher_par%want_re_optical_depth = ini_read_logical( 'param[tau]' ,.false. )
563  fp%fisher_par%want_bias = ini_read_logical( 'param[bias]' ,.false. )
564  fp%fisher_par%want_alpha_SN = ini_read_logical( 'param[alpha_SN]',.false. )
565  fp%fisher_par%want_beta_SN = ini_read_logical( 'param[beta_SN]' ,.false. )
566  fp%fisher_par%want_M0_SN = ini_read_logical( 'param[M0_SN]' ,.false. )
567 
568 #ifdef COSMICFISH_CAMB
569  fp%fisher_par%want_w0_ppf = ini_read_logical( 'param[w0_ppf]' ,.false. )
570  fp%fisher_par%want_wa_ppf = ini_read_logical( 'param[wa_ppf]' ,.false. )
571  fp%fisher_par%want_cs_ppf = ini_read_logical( 'param[cs_ppf]' ,.false. )
572 #endif
573 
574  ! read Cls parameters:
575  if ( fp%cosmicfish_want_cls ) then
576 
577  fp%fisher_cls%Fisher_want_CMB_T = ini_read_logical( 'Fisher_want_CMB_T' , .false. )
578  fp%fisher_cls%Fisher_want_CMB_E = ini_read_logical( 'Fisher_want_CMB_E' , .false. )
579  fp%fisher_cls%Fisher_want_CMB_B = ini_read_logical( 'Fisher_want_CMB_B' , .false. )
580  fp%fisher_cls%Fisher_want_CMB_lensing = ini_read_logical( 'Fisher_want_CMB_lensing' , .false. )
581  fp%fisher_cls%Fisher_want_LSS_lensing = ini_read_logical( 'Fisher_want_LSS_lensing' , .false. )
582  fp%fisher_cls%Fisher_want_LSS_counts = ini_read_logical( 'Fisher_want_LSS_counts' , .false. )
583  fp%fisher_cls%Fisher_want_XC = ini_read_logical( 'Fisher_want_XC' , .false. )
584 
585  fp%fisher_cls%CMB_n_channels = ini_read_int( 'CMB_n_channels' , 0 )
586  fp%fisher_cls%CMB_TT_fsky = ini_read_double( 'CMB_TT_fsky' , 1._dl )
587  fp%fisher_cls%CMB_EE_fsky = ini_read_double( 'CMB_EE_fsky' , 1._dl )
588  fp%fisher_cls%CMB_BB_fsky = ini_read_double( 'CMB_BB_fsky' , 1._dl )
589 
590  fp%fisher_cls%l_max_TT = ini_read_int( 'l_max_TT' , 0 )
591  fp%fisher_cls%l_max_EE = ini_read_int( 'l_max_EE' , 0 )
592  fp%fisher_cls%l_max_BB = ini_read_int( 'l_max_BB' , 0 )
593 
594  allocate( fp%fisher_cls%CMB_temp_sens( fp%fisher_cls%CMB_n_channels ), &
595  & fp%fisher_cls%CMB_pol_sens( fp%fisher_cls%CMB_n_channels ), &
596  & fp%fisher_cls%CMB_fwhm( fp%fisher_cls%CMB_n_channels ) )
597 
598  do i=1, fp%fisher_cls%CMB_n_channels
599  write (numstr,*) i
600  numstr=adjustl(numstr)
601  fp%fisher_cls%CMB_temp_sens(i) = ini_read_double( 'CMB_temp_sens('//trim(numstr)//')', 0._dl )
602  fp%fisher_cls%CMB_pol_sens(i) = ini_read_double( 'CMB_pol_sens('//trim(numstr)//')', 0._dl )
603  fp%fisher_cls%CMB_fwhm(i) = ini_read_double( 'CMB_fwhm('//trim(numstr)//')', 0._dl )
604  end do
605 
606  fp%fisher_cls%LSS_number_windows = ini_read_int( 'LSS_number_windows' , 0 )
607 
608  allocate( fp%fisher_cls%LSS_num_galaxies( fp%fisher_cls%LSS_number_windows ), &
609  & fp%fisher_cls%LSS_intrinsic_ellipticity( fp%fisher_cls%LSS_number_windows ), &
610  & fp%fisher_cls%LSS_fsky( fp%fisher_cls%LSS_number_windows ), &
611  & fp%fisher_cls%LSS_lmax( fp%fisher_cls%LSS_number_windows ) )
612 
613  do i=1, fp%fisher_cls%LSS_number_windows
614  write (numstr,*) i
615  numstr=adjustl(numstr)
616  fp%fisher_cls%LSS_num_galaxies(i) = ini_read_double('LSS_num_galaxies('//trim(numstr)//')', 0._dl )
617  fp%fisher_cls%LSS_intrinsic_ellipticity(i) = ini_read_double('LSS_intrinsic_ellipticity('//trim(numstr)//')', 0._dl )
618  fp%fisher_cls%LSS_fsky(i) = ini_read_double('LSS_fsky('//trim(numstr)//')', 0._dl )
619  fp%fisher_cls%LSS_lmax(i) = ini_read_double('LSS_lmax('//trim(numstr)//')', 0._dl )
620  end do
621 
622  fp%fisher_cls%window_type = ini_read_int( 'window_type' , 1 )
623  fp%fisher_cls%window_alpha = ini_read_double('window_alpha' , 0._dl )
624  fp%fisher_cls%window_beta = ini_read_double('window_beta' , 0._dl )
625  fp%fisher_cls%redshift_zero = ini_read_double('redshift_zero', 0._dl )
626  fp%fisher_cls%photoz_error = ini_read_double('photoz_error' , 0._dl )
627 
628  if ( fp%fisher_par%want_bias ) then
629  ! number of GC windows:
630  ind = 0
631  do i = 1, fp%fisher_cls%LSS_number_windows
632  if ( redshift_w(i)%kind == window_counts ) then
633  ind = ind +1
634  end if
635  end do
636  ! allocate the bias array:
637  allocate( fp%fisher_cls%bias( ind ) )
638  ! write the fiducial:
639  ind = 0
640  do i = 1, fp%fisher_cls%LSS_number_windows
641  if ( redshift_w(i)%kind == window_counts ) then
642  ind = ind +1
643  fp%fisher_cls%bias( ind ) = redshift_w(i)%bias
644  end if
645  end do
646  end if
647 
648  end if
649 
650  ! read SN parameters:
651 
652  ! fiducial SN parameters:
653  fp%fisher_SN%alpha_SN = ini_read_double('alpha_SN' , 0._dl )
654  fp%fisher_SN%beta_SN = ini_read_double('beta_SN' , 0._dl )
655  fp%fisher_SN%M0_SN = ini_read_double('M0_SN' , 0._dl )
656 
657  if ( fp%cosmicfish_want_SN ) then
658 
659  fp%fisher_SN%number_SN_windows = ini_read_int( 'number_SN_windows', 0 )
660 
661  allocate( fp%fisher_SN%SN_number(fp%fisher_SN%number_SN_windows), &
662  & fp%fisher_SN%SN_redshift_start(fp%fisher_SN%number_SN_windows), &
663  & fp%fisher_SN%SN_redshift_end(fp%fisher_SN%number_SN_windows) )
664 
665  do i=1, fp%fisher_SN%number_SN_windows
666  write (numstr,*) i
667  numstr=adjustl(numstr)
668  fp%fisher_SN%SN_number(i) = ini_read_int( 'SN_number('//trim(numstr)//')', 0 )
669  fp%fisher_SN%SN_redshift_start(i) = ini_read_double('SN_redshift_start('//trim(numstr)//')', 0._dl )
670  fp%fisher_SN%SN_redshift_end(i) = ini_read_double('SN_redshift_end('//trim(numstr)//')', 0._dl )
671  end do
672 
673  fp%fisher_SN%total_SN_number = sum( fp%fisher_SN%SN_number )
674 
675  ! SN Fisher Monte Carlo average:
676  fp%fisher_SN%SN_Fisher_MC_samples = ini_read_int( 'SN_Fisher_MC_samples', 1 )
677 
678  ! parameters to generate the SN mock data:
679  fp%fisher_SN%color_dispersion = ini_read_double('color_dispersion' , 0._dl )
680  fp%fisher_SN%stretch_dispersion = ini_read_double('stretch_dispersion' , 0._dl )
681 
682  ! parameters to generate the mock SN covariance:
683  fp%fisher_SN%magnitude_sigma = ini_read_double('magnitude_sigma' , 0._dl )
684  fp%fisher_SN%c_sigmaz = ini_read_double('c_sigmaz' , 0._dl )
685  fp%fisher_SN%sigma_lens_0 = ini_read_double('sigma_lens_0' , 0._dl )
686 
687  fp%fisher_SN%dcolor_offset = ini_read_double('dcolor_offset' , 0._dl )
688  fp%fisher_SN%dcolor_zcorr = ini_read_double('dcolor_zcorr' , 0._dl )
689  fp%fisher_SN%dshape_offset = ini_read_double('dshape_offset' , 0._dl )
690  fp%fisher_SN%dshape_zcorr = ini_read_double('dshape_zcorr' , 0._dl )
691 
692  fp%fisher_SN%cov_ms_offset = ini_read_double('cov_ms_offset' , 0._dl )
693  fp%fisher_SN%cov_ms_zcorr = ini_read_double('cov_ms_zcorr' , 0._dl )
694  fp%fisher_SN%cov_mc_offset = ini_read_double('cov_mc_offset' , 0._dl )
695  fp%fisher_SN%cov_mc_zcorr = ini_read_double('cov_mc_zcorr' , 0._dl )
696  fp%fisher_SN%cov_sc_offset = ini_read_double('cov_sc_offset' , 0._dl )
697  fp%fisher_SN%cov_sc_zcorr = ini_read_double('cov_sc_zcorr' , 0._dl )
698 
699  end if
700 
701  ! read RD parameters:
702  if ( fp%cosmicfish_want_RD ) then
703  fp%fisher_RD%exptype = ini_read_int( 'RD_exptype', 0 )
704  if (fp%fisher_RD%exptype.eq.1) then
705  fp%fisher_RD%number_RD_redshifts = ini_read_int( 'number_RD_redshifts', 0 )
706  allocate(fp%fisher_RD%RD_redshift(fp%fisher_RD%number_RD_redshifts),fp%fisher_RD%RD_number(fp%fisher_RD%number_RD_redshifts))
707  do i=1,fp%fisher_RD%number_RD_redshifts
708  write(red_ind,*) i
709  fp%fisher_RD%RD_redshift(i) = ini_read_double( 'RD_redshift('//trim(adjustl(red_ind))//')',0._dl)
710  fp%fisher_RD%RD_number(i) = ini_read_int( 'RD_source_number('//trim(adjustl(red_ind))//')', 0 ) !MM can be improved
711  end do
712  fp%fisher_RD%obs_time = ini_read_int( 'delta_time', 0 )
713  fp%fisher_RD%signoise = ini_read_int( 'RD_sig_to_noise', 0 )
714  else if (fp%fisher_RD%exptype.eq.2) then
715  write(0,*)'not implemented yet'
716  stop
717  else if (fp%fisher_RD%exptype.eq.3) then
718  write(0,*)'not implemented yet'
719  stop
720  end if
721  end if
722 
723  ! read derived parameters:
724  if ( fp%cosmicfish_want_derived ) then
725 
726  fp%fisher_der%want_omegab = ini_read_logical( 'param[omegab]' ,.false. )
727  fp%fisher_der%want_omegac = ini_read_logical( 'param[omegac]' ,.false. )
728  fp%fisher_der%want_omegan = ini_read_logical( 'param[omegan]' ,.false. )
729  fp%fisher_der%want_omegav = ini_read_logical( 'param[omegav]' ,.false. )
730  fp%fisher_der%want_omegak = ini_read_logical( 'param[omegak]' ,.false. )
731  fp%fisher_der%want_omegam = ini_read_logical( 'param[omegam]' ,.false. )
732  fp%fisher_der%want_theta = ini_read_logical( 'param[theta]' ,.false. )
733  fp%fisher_der%want_mnu = ini_read_logical( 'param[mnu]' ,.false. )
734  fp%fisher_der%want_zre = ini_read_logical( 'param[zre]' ,.false. )
735  fp%fisher_der%want_neff = ini_read_logical( 'param[neff]' ,.false. )
736 
737  fp%fisher_der%want_sigma8 = ini_read_logical( 'param[sigma8]' ,.false. )
738  fp%fisher_der%want_loghubble = ini_read_logical( 'param[loghubble]' ,.false. )
739  fp%fisher_der%want_logDA = ini_read_logical( 'param[logDA]' ,.false. )
740 
741  fp%fisher_der%FD_num_redshift = ini_read_int( 'FD_num_redshift', 0 )
742 
743  allocate( fp%fisher_der%FD_redshift( fp%fisher_der%FD_num_redshift ) )
744 
745  do i=1, fp%fisher_der%FD_num_redshift
746  write (numstr,*) i
747  numstr=adjustl(numstr)
748 
749  fp%fisher_der%FD_redshift(i) = ini_read_double( 'FD_redshift('//trim(numstr)//')', 0._dl )
750  end do
751 
752  end if
753 
754  ! read Mpk parameters:
755  if ( fp%cosmicfish_want_Mpk ) then
756 
757  end if
758 
759  ! dump read parameters to file:
760  if ( present(param_out_name) ) then
761  call ini_savereadvalues(trim(param_out_name),1)
762  else if ( fp%outroot /= '' ) then
763  if (filename /= trim(fp%outroot) //'params.ini') then
764  call ini_savereadvalues(trim(fp%outroot) //'params.ini',1)
765  else
766  write(*,*) 'WARNING: Output _params.ini not created as would overwrite input parameters.'
767  end if
768  end if
769 
770  ! finalize:
771  call ini_close
772  if (.not. camb_validateparams(p)) stop 'Stopped due to parameter error'
773  if (global_error_flag/=0) stop 'Error in initialization.'
774 
775  ! check the consistency of CosmicFish parameters
776  call check_parameters_consistency(p, fp)
777 
778  end subroutine init_cosmicfish_from_file
779 
780  ! ---------------------------------------------------------------------------------------------
782  subroutine check_parameters_consistency( P, FP )
783 
784  implicit none
785 
786  Type(cambparams) :: P
787  Type(cosmicfish_params) :: FP
788 
789  ! general check:
790  if ( .not. fp%cosmicfish_want_cls .and. &
791  & .not. fp%cosmicfish_want_SN .and. &
792  & .not. fp%cosmicfish_want_RD .and. &
793  & .not. fp%cosmicfish_want_Mpk .and. &
794  & .not. fp%cosmicfish_want_derived ) then
795  write(*,*) 'WARNING: you want no Fisher!.'
796  end if
797 
798  ! check the cls Fisher matrix:
799  if ( fp%cosmicfish_want_cls ) then
800  ! l_max check:
801  if ( p%Max_l+20 < fp%fisher_cls%l_max_TT .or. &
802  & p%Max_l+20 < fp%fisher_cls%l_max_EE .or. &
803  & p%Max_l+20 < fp%fisher_cls%l_max_BB ) then
804  write(*,*) 'ERROR: l_max+20 is smaller that l_max_TT or l_max_EE or l_max_BB.'
805  stop
806  end if
807  end if
808 
809  write(*,*) 'WARNING: check_parameters_consistency: not fully implemented.'
810 
811  end subroutine check_parameters_consistency
812 
813 end module init_from_file
init_from_file::init_cosmicfish_from_file
subroutine, public init_cosmicfish_from_file(P, FP, filename, param_out_name)
This subroutine initializes camb parameters from a parameter file.
Definition: 005_init_from_file.f90:44
init_from_file
This module contains the subroutine and functions to initialize camb and cosmicfish from a parameter ...
Definition: 005_init_from_file.f90:26
cosmicfish_types
This module contains the definitions of derived data types used to store the informations about the f...
Definition: 004_cosmicfish_types.f90:26
init_from_file::check_parameters_consistency
subroutine check_parameters_consistency(P, FP)
This subroutine checks the consistency of the cosmicfish parameters with the camb parameters...
Definition: 005_init_from_file.f90:783
cosmicfish_types::cosmicfish_params
This derived data type contains all the informations that the cosmicfish code needs to run...
Definition: 004_cosmicfish_types.f90:215
Generated on Tue Jun 21 2016 03:18:43 for CosmicFish by   doxygen 1.8.9.1