from numpy import array, log10
import pylab
from pymses.analysis.visualization import *
from pymses import RamsesOutput
from pymses.utils import constants as C

# Ramses data
ioutput = 193
ro = RamsesOutput("/data/Aquarius/output/", ioutput)
parts = ro.particle_source(["mass", "level"])

# Map operator : mass
scal_func = ScalarOperator(lambda dset: dset["mass"])

# Map region
center = [ 0.567811, 0.586055, 0.559156 ]
axes = {"los": array([ -0.172935, 0.977948, -0.117099 ])}

# Map processing
mp = fft_projection.MapFFTProcessor(parts, ro.info)
for axname, axis in axes.items():
	cam  = Camera(center=center, line_of_sight_axis=axis, up_vector="z", region_size=[5.0E-1, 4.5E-1], \
			distance=2.0E-1, far_cut_depth=2.0E-1, map_max_size=512)
	map = mp.process(scal_func, cam, surf_qty=True)
	factor = (ro.info["unit_mass"]/ro.info["unit_length"]**2).express(C.Msun/C.kpc**2)
	scale = ro.info["unit_length"].express(C.Mpc)

#	pylab.imshow(map)
#	pylab.show()
	# Save map into HDF5 file
	mapname = "DM_Sigma_%s_%5.5i"%(axname, ioutput)
	h5fname = save_map_HDF5(map, cam, map_name=mapname)

	# Plot map into Matplotlib figure/PIL Image
	fig = save_HDF5_to_plot(h5fname, map_unit=("$M_{\odot}.kpc^{-2}$",factor), axis_unit=("Mpc", scale), cmap="Blues", fraction=0.1)
#	pil_img = save_HDF5_to_img(h5fname, cmap="Blues", fraction=0.1)
	
	# Save map into PNG image file
#	save_HDF5_to_plot(h5fname, map_unit=("$M_{\odot}.kpc^{-2}$",factor), \
#			axis_unit=("Mpc", scale), img_path="./", cmap="Blues", fraction=0.1)
#	save_HDF5_to_img(h5fname, img_path="./", cmap="Blues", fraction=0.1)

# pylab.show()