Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
On Thursday, 7th July from 1 to 3 pm there will be a maintenance with a short downtime of GitLab.
Open sidebar
ift
IMAGINE
Commits
7778ae29
Commit
7778ae29
authored
Oct 27, 2017
by
Theo Steininger
Browse files
Adapted IMAGINE to newest version of HamX.
parent
1a2cb14e
Changes
15
Hide whitespace changes
Inline
Side-by-side
docker/Dockerfile
View file @
7778ae29
...
...
@@ -41,14 +41,9 @@ ENV HEALPIX /home/Downloads/Healpix_3.31
RUN
pip
install
healpy
#Hammurabi
# RUN wget https://sourceforge.net/code-snapshots/svn/h/ha/hammurabicode/code/hammurabicode-code-41-trunk.zip
# RUN unzip hammurabicode-code-41-trunk.zip
# WORKDIR hammurabicode-code-41-trunk
RUN
wget https://sourceforge.net/code-snapshots/svn/h/ha/hammurabicode/code/hammurabicode-code-60-trunk.zip
RUN
unzip hammurabicode-code-60-trunk.zip
WORKDIR
hammurabicode-code-60-trunk
COPY
hammurabi_Makefile Makefile
RUN
[
-r
/root/.healpix/3_31_Linux/config
]
&&
.
/root/.healpix/3_31_Linux/config
&&
make hammurabi
&&
make clean
RUN
git clone https://bitbucket.org/ricphy/hamx
-b
generic_makefile
WORKDIR
hamx
RUN
make
WORKDIR
..
#(Py)MultiNest
...
...
@@ -101,5 +96,12 @@ WORKDIR ..
#hampy
RUN
pip
install
jupyter pandas
#IMAGINE
RUN
git clone https://gitlab.mpcdf.mpg.de/ift/IMAGINE.git
-b
master
WORKDIR
IMAGINE
RUN
python setup.py
install
WORKDIR
..
docker/hammurabi_Makefile
deleted
100644 → 0
View file @
1a2cb14e
#--------------------------------------------------------
#--------------------------------------------------------
#
# This is the hammurabi Makefile. The first part needs to be configured for your system.
#
# You can choose to build:
#
# hammurabi: by default, the basic hammurabi, pure C++, for constant or grid TEs and simple analytic models for CREs
# hammurabi.ne2001: including linking to the (smooth) NE2001 Fortran code using cfortran.h
# hammurabi.galprop: including both NE2001 and Galprop to simulate the CRE propagation with the hammurabi magnetic field first.
# hammurabi.debug: including all code, debug flags.
# all: builds all of these.
#
#--------------------------------------------------------
#--------------------------------------------------------
#
# Give the locations of the various libraries. Should containt include and lib subdirectories.
#
#BASE_DIR = $(HOME)/space/sw/
#HAMMURABI_HOME = $(BASE_DIR)/hammurabi/branches/mine
#GSL = $(BASE_DIR)/gsl/gsl-1.13/build/
#FFTW = $(BASE_DIR)/fftw/fftw-3.2.2/build/
#CFITSIO = $(BASE_DIR)/cfitsio/cfitsio/build/
BASE_DIR = /home/Downloads/hammurabicode-code
HAMMURABI_HOME = $(BASE_DIR)
GSL = /usr/local
FFTW = /usr/local
CFITSIO = /usr/local
#
# These are needed if you build with Galprop
#
GALPROP = $(BASE_DIR)/galprop/galprop/build
CCFITS = $(BASE_DIR)/ccfits/CCfits//build/
CLHEP = $(BASE_DIR)/clhep/clhep/build/
#
# What c++ compiler you are using?
#
CXX = g++
#
# On OSX, maybe
#CXX = g++-mp-5
#
# Compile options? 64bit system (-m64), optimization (-O2), openMP (-fopenmp), warnings (-Wall), etc.
#
CXXFLAGS = -m64 -O2 -fopenmp -Wall
#
# With OSX g++-mp-4.3
#CXXFLAGS = -fopenmp -O2 -g -fno-inline-functions -Wall -Wextra -Wno-unknown-pragmas -ansi
#
# Which Fortran compiler you are using? Only if compiling with NE2001 and/or Galprop
#
FC = gfortran
#
# On OSX, maybe
#FC = gfortran-mp-5
# For gfortran
FCFLAGS = -O2 -ffixed-form -ffixed-line-length-132 -c
#
# For ifort
#FCFLAGS =-O -extend-source -c
#
# For f77 ?
#FCFLAGS =
#
# For linking to C++ with gcc
#
LFORTFLAGS= -Df2cFortran
#
# Additional options depending on how you want to run:
#
# Use openMP. Right now, this is mandatory. (Contact
# trjaffe@gmail.com for help compiling without if for some reason you
# need it.)
LFLAGS_CUSTOM = -fopenmp
#
# OSX with gfortran-mp-4.3
#LFLAGS_CUSTOM = -fopenmp -lgcc_s.1
#
#
#
# For J. West using f77 and gcc 4.1.2?
#
#LFORTFLAGS += -L/usr/lib/gcc/x86_64-redhat-linux/3.4.6/ -lgcc -lg2c
#LFLAGS_CUSTOM += -lgcc_s
#------------------------------------------------------------
#------------------------------------------------------------
#
# From here onwards, nothing should need to be changed.
#
#------------------------------------------------------------
#------------------------------------------------------------
hammurabi.debug: CXXFLAGS=-m64 -g -Wall
hammurabi: ALL_INC = -I$(CFITSIO)/include -I$(GSL)/include -I. -I$(FFTW)/include -I$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/include/
hammurabi.ne2001: ALL_INC = -I$(CFITSIO)/include -I$(GSL)/include -I. -I$(FFTW)/include -I$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/include/
hammurabi.galprop: ALL_INC = -I$(CFITSIO)/include -I$(GALPROP)/.. -I$(GALPROP)/include/ -I$(CCFITS)/include/ -I$(CLHEP)/include -I$(GSL)/include -I. -I$(FFTW)/include -I$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/include/
hammurabi.debug: ALL_INC = -I$(CFITSIO)/include -I$(GALPROP)/.. -I$(GALPROP)/include/ -I$(CCFITS)/include/ -I$(CLHEP)/include -I$(GSL)/include -I. -I$(FFTW)/include -I$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/include/
libhamprop.a: ALL_INC = -I$(CFITSIO)/include -I$(GALPROP)/.. -I$(GALPROP)/include/ -I$(GSL)/include -I. -I$(FFTW)/include -I$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/include/
hammurabi: ALL_L = -L$(CFITSIO)/lib -L$(GSL)/lib -L$(FFTW)/lib -L$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/lib/ -L.
hammurabi.ne2001: ALL_L = -L$(CFITSIO)/lib -L$(GSL)/lib -L$(FFTW)/lib -L$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/lib/ -L.
hammurabi.galprop: ALL_L = -L$(CFITSIO)/lib -L$(GALPROP)/lib/ -L$(CCFITS)/lib/ -L$(CLHEP)/lib -L$(GSL)/lib -L$(FFTW)/lib -L$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/lib/ -L.
hammurabi.debug: ALL_L = -L$(CFITSIO)/lib -L$(GALPROP)/lib/ -L$(CCFITS)/lib/ -L$(CLHEP)/lib -L$(GSL)/lib -L$(FFTW)/lib -L$(HEALPIX)/src/cxx/$(HEALPIX_TARGET)/lib/ -L.
hammurabi : CXXFLAGS+= $(ALL_INC) -c
hammurabi.ne2001: CXXFLAGS+= $(ALL_INC) -c -DNE2001
hammurabi.galprop: CXXFLAGS+= $(ALL_INC) -c -DGALDEF_PATH=\"./GALDEF\" -DFITSDATA_PATH=\"./FITS/\" -DDATA_PATH=\"./DATA\" -DGALPROP -DNE2001
hammurabi.debug : CXXFLAGS+= $(ALL_INC) -c -DGALDEF_PATH=\"./GALDEF\" -DFITSDATA_PATH=\"./FITS/\" -DDATA_PATH=\"./DATA\" -DGALPROP -DNE2001 -g -Wall
# Turns off sanity check to avoid some circularity. See Galprop README.tess
libhamprop.a: CXXFLAGS+= $(ALL_INC) -c -DGALDEF_PATH=\"./GALDEF\" -DFITSDATA_PATH=\"./FITS/\" -DDATA_PATH=\"./DATA\" -DLIBHAMPROP -DGALPROP -DNE2001
hammurabi: LFLAGS = $(ALL_L) -lhammurabi -lhealpix_cxx -lcxxsupport -lc_utils -lfftpack -lcfitsio -lgsl -lgslcblas -lm -lfftw3 -lfftw3_omp $(LFLAGS_CUSTOM)
hammurabi.ne2001: LFLAGS = $(ALL_L) -lhammurabi -lNE2001 -lhealpix_cxx -lcxxsupport -lc_utils -lfftpack -lcfitsio -lgsl -lgslcblas -lm -lfftw3 -lfftw3_omp -lgfortran $(LFLAGS_CUSTOM)
hammurabi.debug: LFLAGS = $(ALL_L) -lhammurabi -lNE2001 -lhealpix_cxx -lcxxsupport -lc_utils -lfftpack -lcfitsio -lgsl -lgslcblas -lm -lfftw3 -lfftw3_omp -lgfortran -lgalprop -lskymap -lCCfits -lCLHEP $(LFLAGS_CUSTOM)
hammurabi.galprop: LFLAGS = $(ALL_L) -lhammurabi -lNE2001 -lhealpix_cxx -lcxxsupport -lc_utils -lfftpack -lcfitsio -lgsl -lgslcblas -lm -lfftw3 -lfftw3_omp -lgfortran -lgalprop -lskymap -lCCfits -lCLHEP $(LFLAGS_CUSTOM)
default:
make hammurabi
all:
make clean ; make hammurabi ; make clean ; make hammurabi.ne2001; make clean; make hammurabi.galprop ; make clean ; make hammurabi.debug
%.o : %.cpp
$(CXX) $(CXXFLAGS) -o $@ $<
%.o : %.cc
$(CXX) $(CXXFLAGS) -o $@ $<
%.o : %.c
$(CXX) $(CXXFLAGS) -o $@ $<
%.o : %.f
$(FC) $(FCFLAGS) -o $@ $<
class_TE_density.o: class_TE_density.cpp
$(CXX) $(CXXFLAGS) $(LFORTFLAGS) -o $@ $<
NE2001_OBJ = dmdsm.NE2001.o density.NE2001.o neclumpN.o nevoidN.o neLISM.NE2001.o
HAMMURABI_OBJ = class_Integrator.o class_List.o class_B_field2.o tess_tools.o class_TE_density.o class_CRE.o namespace_Vec_Handling.o class_Dust.o
run:
mkdir run
$(HAMMURABI_OBJ): hammurabi.h CGS_units_file.h proto*.h tess_tools.h
libNE2001.a: $(NE2001_OBJ)
ar rc $@ $(NE2001_OBJ)
ranlib $@
hammurabi: $(HAMMURABI_OBJ) libhammurabi.a run hammurabi.o
$(CXX) -o run/hammurabi hammurabi.o $(LFLAGS) $(LFORTFLAGS)
hammurabi.ne2001: $(HAMMURABI_OBJ) $(NE2001_OBJ) $(NE2001_DATA) libNE2001.a libhammurabi.a run hammurabi.o
$(CXX) -o run/hammurabi.ne2001 hammurabi.o $(LFLAGS) $(LFORTFLAGS)
hammurabi.galprop: $(HAMMURABI_OBJ) $(NE2001_OBJ) $(NE2001_DATA) libNE2001.a libhammurabi.a run hammurabi.o
$(CXX) -o run/hammurabi.galprop hammurabi.o $(LFLAGS) $(LFORTFLAGS)
hammurabi.debug: $(HAMMURABI_OBJ) $(NE2001_DATA) libNE2001.a libhammurabi.a run hammurabi.o
$(CXX) -o run/hammurabi.debug hammurabi.o $(LFLAGS) $(LFORTFLAGS)
ifeq ($(OSTYPE),darwin)
install_name_tool -change "@rpath/libCLHEP-2.3.3.1.dylib" "$(CLHEP)/lib/libCLHEP-2.3.3.1.dylib" run/hammurabi.debug
endif
print_constants: print_constants.o
$(CXX) -o run/print_constants print_constants.o
# To link in B_field to galprop, use this:
libhamprop.a: class_B_field2.o tess_tools.o namespace_Vec_Handling.o
ar cru libhamprop.a class_B_field2.o tess_tools.o namespace_Vec_Handling.o
libhammurabi.a: $(HAMMURABI_OBJ)
ar cru libhammurabi.a $(HAMMURABI_OBJ)
clean:
rm *.a *.o
test:
(cd unit_test; ./test.csh)
test_big:
(cd unit_test; ./test.csh big)
tarfile:
tar cvzf hammurabi.tgz *cpp *h Makefile *.f README hampy
unittar:
tar cvzf hammurabi_unit_test_inputs.tgz unit_test/inputs unit_test/GALDEF unit_test/FITS unit_test/negrid_n400.bin
tar cvzf hammurabi_unit_test_ref.tgz unit_test/ref unit_test/ref.big
tar cvzf hammurabi_unit_test_ref.mini.tgz unit_test/ref.mini
tar cvzf hammurabi_unit_test_ref.big.tgz unit_test/ref.big
imagine/likelihoods/ensemble_likelihood/ensemble_likelihood.py
View file @
7778ae29
...
...
@@ -50,14 +50,16 @@ class EnsembleLikelihood(Likelihood):
u_val
=
obs_val
-
obs_mean
# compute quantities for OAS estimator
mu
=
np
.
vdot
(
u_val
,
u_val
)
*
weight
/
n
mu
=
np
.
vdot
(
u_val
,
u_val
)
*
weight
/
k
self
.
logger
.
debug
(
"mu: %f"
%
mu
)
alpha
=
(
np
.
einsum
(
u_val
,
[
0
,
1
],
u_val
,
[
2
,
1
])
**
2
).
sum
()
# correct the volume factor: one factor comes from the internal scalar
# product and one from the trace
alpha
*=
weight
**
2
numerator
=
alpha
+
mu
**
2
denominator
=
(
k
+
1
)
*
(
alpha
-
(
mu
**
2
)
/
n
)
numerator
=
(
1
-
2.
/
n
)
*
alpha
+
mu
**
2
denominator
=
(
k
+
1
-
2.
/
n
)
*
(
alpha
-
(
mu
**
2
)
/
n
)
if
denominator
==
0
:
rho
=
1
...
...
@@ -77,7 +79,7 @@ class EnsembleLikelihood(Likelihood):
"DiagonalOperator."
)
A_bare_diagonal
=
data_covariance_operator
.
diagonal
(
bare
=
True
)
A_bare_diagonal
.
val
+=
rho
*
mu
A_bare_diagonal
.
val
+=
rho
*
mu
/
n
A
=
DiagonalOperator
(
domain
=
data_covariance_operator
.
domain
,
diagonal
=
A_bare_diagonal
,
...
...
imagine/magnetic_fields/wmap3yr_magnetic_field/wmap3yr_magnetic_field.py
View file @
7778ae29
...
...
@@ -8,13 +8,13 @@ class WMAP3yrMagneticField(MagneticField):
@
property
def
descriptor_lookup
(
self
):
lookup
=
\
{
'b0'
:
(
'./Galaxy/MagneticField/Regular/WMAP/b0'
,
'value'
)
,
'psi0'
:
(
'./Galaxy/MagneticField/Regular/WMAP/psi0'
,
'value'
)
,
'psi1'
:
(
'./Galaxy/MagneticField/Regular/WMAP/psi1'
,
'value'
)
,
'chi0'
:
(
'./Galaxy/MagneticField/Regular/WMAP/chi0'
,
'value'
)
,
'random_rms'
:
(
'./Galaxy/MagneticField/Random/Iso/rms'
,
'value'
)
,
'random_rho'
:
(
'./Galaxy/MagneticField/Random/Anisoglob/rho'
,
'value'
)
}
{
'b0'
:
[
'./Galaxy/MagneticField/Regular/WMAP/b0'
,
'value'
]
,
'psi0'
:
[
'./Galaxy/MagneticField/Regular/WMAP/psi0'
,
'value'
]
,
'psi1'
:
[
'./Galaxy/MagneticField/Regular/WMAP/psi1'
,
'value'
]
,
'chi0'
:
[
'./Galaxy/MagneticField/Regular/WMAP/chi0'
,
'value'
]
,
'random_rms'
:
[
'./Galaxy/MagneticField/Random/Iso/rms'
,
'value'
]
,
'random_rho'
:
[
'./Galaxy/MagneticField/Random/Anisoglob/rho'
,
'value'
]
}
return
lookup
def
_create_field
(
self
):
...
...
imagine/observers/hammurapy/hammurapy.py
View file @
7778ae29
# -*- coding: utf-8 -*-
import
abc
import
os
import
tempfile
import
subprocess
import
xml.etree.ElementTree
as
et
import
healpy
import
numpy
as
np
from
d2o
import
distributed_data_object
from
nifty
import
HPSpace
from
imagine.observers.observer
import
Observer
from
imagine.magnetic_fields.magnetic_field
import
MagneticField
from
.observable_mixins
import
ObservableMixin
from
.model_mixins
import
MagneticFieldModel
from
imagine.observables
import
Observable
class
Hammurapy
(
Observer
):
def
__init__
(
self
,
hammurabi_executable
,
input_directory
=
'./input'
,
working_directory_base
=
'.'
,
nside
=
64
):
def
__init__
(
self
,
hammurabi_executable
,
magnetic_field_model
,
observables
,
input_directory
=
'./input'
,
working_directory_base
=
'.'
,
nside
=
64
):
self
.
hammurabi_executable
=
os
.
path
.
abspath
(
hammurabi_executable
)
if
not
isinstance
(
magnetic_field_model
,
MagneticFieldModel
):
raise
TypeError
(
"magnetic_field_model must be an instance of the "
"MagneticField class."
)
self
.
magnetic_field_model
=
magnetic_field_model
if
not
isinstance
(
observables
,
list
):
if
isinstance
(
observables
,
tuple
):
observables
=
list
(
observables
)
else
:
observables
=
[
observables
]
for
obs
in
observables
:
if
not
isinstance
(
obs
,
ObservableMixin
):
raise
TypeError
(
"observables must be an instance of the "
"ObservableMixin class."
)
self
.
observables
=
observables
self
.
input_directory
=
os
.
path
.
abspath
(
input_directory
)
self
.
working_directory_base
=
os
.
path
.
abspath
(
working_directory_base
)
self
.
nside
=
int
(
nside
)
self
.
_hpSpace
=
HPSpace
(
nside
=
self
.
nside
)
self
.
last_call_log
=
""
@
abc
.
abstract
property
@
property
def
magnetic_field_class
(
self
):
return
M
agnetic
F
ield
return
self
.
magnetic_field_model
.
m
agnetic
_f
ield
_class
def
_make_temp_folder
(
self
):
prefix
=
os
.
path
.
join
(
self
.
working_directory_base
,
'temp_hammurabi_'
)
...
...
@@ -50,21 +71,6 @@ class Hammurapy(Observer):
else
:
self
.
logger
.
warning
(
'Could not delete %s'
%
path
)
def
_read_fits_file
(
self
,
path
,
name
,
nside
):
map_path
=
os
.
path
.
join
(
path
,
name
)
result_list
=
[]
i
=
0
while
True
:
try
:
loaded_map
=
healpy
.
read_map
(
map_path
,
verbose
=
False
,
field
=
i
)
# loaded_map = healpy.ud_grade(loaded_map, nside_out=nside)
result_list
+=
[
loaded_map
]
i
+=
1
except
IndexError
:
break
return
result_list
def
_call_hammurabi
(
self
,
path
):
temp_process
=
subprocess
.
Popen
(
[
self
.
hammurabi_executable
,
'parameters.xml'
],
...
...
@@ -72,12 +78,31 @@ class Hammurapy(Observer):
cwd
=
path
)
self
.
last_call_log
=
temp_process
.
communicate
()[
0
]
def
_initialize_observable_dict
(
self
,
observable_dict
,
magnetic_field
):
pass
def
_initialize_observable_dict
(
self
,
magnetic_field
):
observable_dict
=
{}
for
observable
in
self
.
observables
:
observable_dict
=
self
.
_initialize_observable_dict_helper
(
observable_dict
,
magnetic_field
,
observable
)
return
observable_dict
def
_initialize_observable_dict_helper
(
self
,
observable_dict
,
magnetic_field
,
observable
):
ensemble_space
=
magnetic_field
.
domain
[
0
]
for
component
in
observable
.
component_names
:
# It is important to initialize the Observables with an explicit
# value. Otherwise the d2o will not instantaneuosly be created
# (c.f. lazy object creation).
observable_dict
[
component
]
=
Observable
(
val
=
0
,
domain
=
(
ensemble_space
,
self
.
_hpSpace
),
distribution_strategy
=
'equal'
)
return
observable_dict
def
_build_parameter_dict
(
self
,
parameter_dict
,
magnetic_field
,
local_ensemble_index
):
parameter_dict
.
update
(
self
.
magnetic_field_model
.
parameter_dict
)
parameter_dict
.
update
(
{(
'./Interface/fe_grid'
,
'read'
):
'1'
,
(
'./Interface/fe_grid'
,
'filename'
):
...
...
@@ -107,7 +132,7 @@ class Hammurapy(Observer):
parameter_dict
.
update
(
{(
'./Grid/Integration/nside'
,
'value'
):
self
.
nside
})
def
_write_parameter_
dict
(
self
,
parameter_dict
,
working_directory
):
def
_write_parameter_
file
(
self
,
working_directory
):
# load the default xml
try
:
default_parameters_xml
=
os
.
path
.
join
(
self
.
input_directory
,
...
...
@@ -129,8 +154,78 @@ class Hammurapy(Observer):
'parameters.xml'
)
tree
.
write
(
parameters_file_path
)
def
_write_parameter_xml
(
self
,
magnetic_field
,
local_ensemble_index
,
working_directory
):
# load the default xml
try
:
default_parameters_xml
=
os
.
path
.
join
(
self
.
input_directory
,
'default_parameters.xml'
)
tree
=
et
.
parse
(
default_parameters_xml
)
except
IOError
:
import
imagine
module_path
=
os
.
path
.
split
(
imagine
.
observers
.
hammurapy
.
__file__
)[
0
]
default_parameters_xml
=
os
.
path
.
join
(
module_path
,
'input/default_parameters.xml'
)
tree
=
et
.
parse
(
default_parameters_xml
)
root
=
tree
.
getroot
()
# modify the default_parameters.xml
custom_parameters
=
[
[
'./Grid/Integration/shell/auto/nside_min'
,
'value'
,
self
.
nside
],
[
'./Interface/fe_grid'
,
'read'
,
'1'
],
[
'./Interface/fe_grid'
,
'filename'
,
os
.
path
.
join
(
self
.
input_directory
,
'fe_grid.bin'
)]
]
# access the magnetic-field's random-eed d2o directly, since we
# know that the distribution strategy is the same for the
# randam samples and the magnetic field itself
random_seed
=
magnetic_field
.
random_seed
.
data
[
local_ensemble_index
]
custom_parameters
+=
[[
'./Galaxy/MagneticField/Random'
,
'seed'
,
random_seed
]]
for
key
,
value
in
magnetic_field
.
parameters
.
iteritems
():
desc
=
magnetic_field
.
descriptor_lookup
[
key
]
custom_parameters
+=
[
desc
+
[
value
]]
# set up grid parameters
grid_space
=
magnetic_field
.
domain
[
1
]
lx
,
ly
,
lz
=
\
np
.
array
(
grid_space
.
shape
)
*
np
.
array
(
grid_space
.
distances
)
/
2.
nx
,
ny
,
nz
=
grid_space
.
shape
custom_parameters
+=
[[
'./Grid/Box/x_min'
,
'value'
,
-
lx
],
[
'./Grid/Box/x_max'
,
'value'
,
lx
],
[
'./Grid/Box/y_min'
,
'value'
,
-
ly
],
[
'./Grid/Box/y_max'
,
'value'
,
ly
],
[
'./Grid/Box/z_min'
,
'value'
,
-
lz
],
[
'./Grid/Box/z_max'
,
'value'
,
lz
],
[
'./Grid/Box/nx'
,
'value'
,
nx
],
[
'./Grid/Box/ny'
,
'value'
,
ny
],
[
'./Grid/Box/nz'
,
'value'
,
nz
]]
for
parameter
in
custom_parameters
:
root
.
find
(
parameter
[
0
]).
set
(
parameter
[
1
],
str
(
parameter
[
2
]))
self
.
magnetic_field_model
.
update_parameter_xml
(
root
)
for
observable
in
self
.
observables
:
observable
.
update_parameter_xml
(
root
)
parameters_file_path
=
os
.
path
.
join
(
working_directory
,
'parameters.xml'
)
tree
.
write
(
parameters_file_path
)
def
_fill_observable_dict
(
self
,
observable_dict
,
working_directory
,
ensemble_index
):
local_ensemble_index
):
for
observable
in
self
.
observables
:
observable
.
fill_observable_dict
(
observable_dict
=
observable_dict
,
working_directory
=
working_directory
,
local_ensemble_index
=
local_ensemble_index
,
nside
=
self
.
nside
)
return
observable_dict
def
__call__
(
self
,
magnetic_field
):
...
...
@@ -139,9 +234,8 @@ class Hammurapy(Observer):
raise
ValueError
(
"Given magnetic field is not a subclass of"
+
" %s"
%
str
(
self
.
magnetic_field_class
))
observable_dict
=
{}
self
.
_initialize_observable_dict
(
observable_dict
=
observable_dict
,
magnetic_field
=
magnetic_field
)
observable_dict
=
self
.
_initialize_observable_dict
(
magnetic_field
=
magnetic_field
)
# iterate over ensemble and put result into result_observable
# get the local shape by creating a dummy d2o
...
...
@@ -150,8 +244,6 @@ class Hammurapy(Observer):
distribution_strategy
=
'equal'
,
dtype
=
np
.
float
)
parameter_dict
=
{}
local_length
=
dummy
.
distributor
.
local_length
for
local_ensemble_index
in
xrange
(
local_length
):
self
.
logger
.
debug
(
"Processing local_ensemble_index %i."
%
...
...
@@ -159,13 +251,10 @@ class Hammurapy(Observer):
# create a temporary folder
working_directory
=
self
.
_make_temp_folder
()
self
.
_build_parameter_dict
(
parameter_dict
=
parameter_dict
,
self
.
_write_parameter_xml
(
magnetic_field
=
magnetic_field
,
local_ensemble_index
=
local_ensemble_index
)
self
.
_write_parameter_dict
(
parameter_dict
=
parameter_dict
,
working_directory
=
working_directory
)
local_ensemble_index
=
local_ensemble_index
,
working_directory
=
working_directory
)
# call hammurabi
self
.
_call_hammurabi
(
working_directory
)
...
...
@@ -180,6 +269,7 @@ class Hammurapy(Observer):
self
.
last_call_log
)
raise
finally
:
self
.
_remove_folder
(
working_directory
)
pass
#self._remove_folder(working_directory)
return
observable_dict
imagine/observers/hammurapy/input/default_parameters.xml
View file @
7778ae29
<?xml version="1.0"?>
<!-- FULL PARAMETER-SET FOR ALL MODULES -->
<!-- "cue" acts as switches for usage -->
<!-- author, Jiaxin Wang -->
<!-- email, jiwang@sissa.it -->
<!-- "cue", switches for usage -->
<root>
<!-- final output .fits maps -->
<Output>
<!-- Dispersion Measure -->
<DM
cue=
"0"
filename=
"dm.fits"
/>
<!-- Synchrotron Emission -->
<Sync
cue=
"0"
filename=
"sync.fits"
/>
<!-- Faraday Depth -->
<Faraday
cue=
"0"
filename=
"fd.fits"
/>
<!-- <DM cue="1" filename="dm.fits"/> -->
<!-- synchrotron dmission -->
<!-- frequency @ GHz -->
<!-- <Sync cue="1" freq="23" filename="iqu_sync_23.fits"/> -->
<!-- <Sync cue="1" freq="2.4" filename="iqu_sync_2.4.fits"/> -->
<!-- Faraday depth -->
<!-- <Faraday cue="1" filename="fd.fits"/> -->
</Output>
<!-- INTERACTIVE IN/OUTPUT -->
<!-- physical field grid in/out -->
<!-- resolution defined in './Grid/Box' -->
<!-- CRE field resolution defined in './CRE/Numeric' -->
<Interface>
<!-- RESOLUTION DEFINED IN './Grid/Box' -->
<!-- regular magnetic field -->
<breg_grid
read=
"0"
write=
"0"
filename=
"breg.bin"
/>
<!-- turbulent magnetic field -->
...
...
@@ -23,13 +29,12 @@
<fe_grid
read=
"0"
write=
"0"
filename=
"ymw16.bin"
/>
<!-- turbulent free electron field -->
<fernd_grid
read=
"0"
write=
"0"
filename=
"fernd.bin"
/>
<!-- cosmic ray electron -->
<!-- cosmic ray electron
field
-->
<cre_grid
read=
"0"
write=
"0"
filename=
"cre.bin"
/>
</Interface>
<!-- HAMMURABI GRID -->
<Grid
type=
"3D"
ec_frame=
"0"
>
<!-- ec_frame offers a local simulation timesaving scheme -->
<!-- main field grid -->
<Grid
type=
"3D"
>
<SunPosition>
<x
value=
"-8.3"
/>
<!-- kpc -->
<y
value=
"0"
/>
<!-- kpc -->
...
...
@@ -37,34 +42,49 @@
</SunPosition>
<Box>
<!-- grid vertex number
in each direction
-->
<!-- grid vertex number -->
<nx
value=
"400"
/>
<ny
value=
"400"
/>
<nz
value=
"160"
/>
<!-- box length in each direction -->
<lx
value=
"40.0"
/>
<!-- kpc -->
<ly
value=
"40.0"
/>
<!-- kpc -->
<lz
value=
"8.0"
/>
<!-- kpc -->
<nz
value=
"80"
/>
<!-- box limit, in Galactic-centric frame -->
<x_min
value=
"-20.0"
/>
<!-- kpc -->
<x_max
value=
"20.0"
/>
<!-- kpc -->
<y_min
value=
"-20.0"
/>
<!-- kpc -->
<y_max
value=
"20.0"
/>
<!-- kpc -->
<z_min
value=
"-4.0"
/>
<!-- kpc -->
<z_max
value=
"4.0"
/>
<!-- kpc -->