from pylab import *
from matplotlib import pyplot


def pHG_analytical(mu,g):
	phg=(1.0-g*g)/(pow((1.0+g*g-2.0*g*mu),3.0/2.0))
	return phg


# parameter
nrofmoments = int(sys.argv[1]) # uebergabeparameter 1
print('l=' + str(nrofmoments))

g = float(sys.argv[2]) # uebergabeparameter 2
print('g=' + str(g))

# save or show? ### dont use show if run by shell-script!!!
save=1 # to show set 0, to save set 1

filename='tmp/phg_g'+str(g)+'_moments'+str(nrofmoments)+'_'+'.dat'
img_filename='tmp/phg_g'+str(g)+'_moments'+str(nrofmoments)+'_'+'.png'

print 'plot ...'

data=loadtxt(filename)
mu=data[:,0]
phg=data[:,1]

pyplot.plot(mu,phg,label='libradtran')
pyplot.plot(mu,pHG_analytical(mu,g),label='analytical')

pyplot.title('g='+str(g)+'    # moments='+str(nrofmoments))
pyplot.xlabel('mu=cos(theta)')
pyplot.ylabel('pHG')

if save==0:
   pyplot.show()

if save==1:
   print 'save plot ...'
   pyplot.savefig(img_filename)



exit()