#!/work/program/EMAN2.1-release/Python/bin/python


import numpy
import string
import emx
from math import *

import sys

def clean_list(param_list, dirty = ''):

	clean_param_list = []	

	for i in range(len(param_list)):
		param = param_list[i]
		if (param == dirty):
			continue
		clean_param_list.append(param)

	return clean_param_list

def read_param(param_set,apix):

	p_num = string.atoi(param_set[0])
        defocus = string.atof(param_set[1])
        xshift = string.atof(param_set[2])
        yshift = string.atof(param_set[3])
        alpha = string.atof(param_set[4])*pi/180.0
        beta = string.atof(param_set[5])*pi/180.0
        gamma = string.atof(param_set[6])*pi/180.0

	t14 = xshift*apix
	t24 = yshift*apix
	t34 = 0.0

#	mat_alpha = [[cos(alpha),-sin(alpha),0],[sin(alpha),cos(alpha),0],[0,0,1]]
	mat_alpha = [[cos(alpha),sin(alpha),0],[-sin(alpha),cos(alpha),0],[0,0,1]]

#	mat_beta = [[1,0,0],[0,-sin(beta),cos(beta)],[0,cos(beta),sin(beta)]]
	mat_beta = [[1,0,0],[0,sin(beta),cos(beta)],[0,cos(beta),-sin(beta)]]

#	mat_gamma = [[cos(gamma),-sin(gamma),0],[sin(gamma),cos(gamma),0],[0,0,1]]
	mat_gamma = [[cos(gamma),sin(gamma),0],[-sin(gamma),cos(gamma),0],[0,0,1]]
	
	[[t11, t12, t13],[t21,t22,t23],[t31,t32,t33]]=numpy.dot(mat_gamma,numpy.dot(mat_beta,mat_alpha))

	return [p_num, defocus, xshift, yshift, alpha, beta, gamma, t11, t12, t13, t14, t21, t22, t23, t24, t31, t32, t33, t34]


def main():

	root_name = 'micrograph'
	emx_file = root_name + '_param.emx'

	apix = 1.42
	acVolt = 300
	boxsize = 200
	micrograph_n = 50

	emxData = emx.EmxData()
	emxData_m = emx.EmxData()
	emxData_p = emx.EmxData()
	xmlMapper = emx.EmxXmlMapper(emxData)

	for i in range(micrograph_n):

		mic_name = root_name + '_' + ('%3d'%(i+1)).replace(' ','0')
		mrcs_param_file = mic_name + '.params'
		mrcs_file = mic_name + '.mrcs'
		mrc_file = mic_name + '.mrc'

		param_list = file(mrcs_param_file)
	      	param = param_list.readline()
		param = param_list.readline()
		param_set = clean_list(param.split(' '))
	
		[p_num, defocus, xshift, yshift, alpha, beta, gamma, t11, t12, t13, t14, t21, t22, t23, t24, t31, t32, t33, t34] = read_param(param_set, apix)
		
		print p_num, defocus, xshift, yshift, alpha, beta, gamma

		m1 = emx.EmxMicrograph(mrc_file)
		m1.set('acceleratingVoltage', acVolt)
		m1.set('defocusU', defocus)
		m1.set('defocusV', defocus)
		m1.set('pixelSpacing__X', apix)
		m1.set('pixelSpacing__Y', apix)
		emxData_m.addObject(m1)
	
		j = 1
		while (len(param) > 0):
			
			if (j > 1):
				param_set = clean_list(param.split(' '))

				[p_num, defocus, xshift, yshift, alpha, beta, gamma, t11, t12, t13, t14, t21, t22, t23, t24, t31, t32, t33, t34] = read_param(param_set, apix)

				print p_num, defocus, xshift, yshift, alpha, beta, gamma

			p1 = emx.EmxParticle(mrcs_file,p_num)
			p1.set('boxSize__X', boxsize)
			p1.set('boxSize__Y', boxsize)
			p1.set('defocusU', defocus)
			p1.set('defocusV', defocus)
			p1.set('defocusUAngle', 0)
			p1.set('pixelSpacing__X', apix)
			p1.set('pixelSpacing__Y', apix)
			p1.set('transformationMatrix__t11',t11)
			p1.set('transformationMatrix__t12',t12)
                        p1.set('transformationMatrix__t13',t13)
                        p1.set('transformationMatrix__t14',t14)
                        p1.set('transformationMatrix__t21',t21)
                        p1.set('transformationMatrix__t22',t22)
                        p1.set('transformationMatrix__t23',t23)
                        p1.set('transformationMatrix__t24',t24)
                        p1.set('transformationMatrix__t31',t31)
                        p1.set('transformationMatrix__t32',t32)
                        p1.set('transformationMatrix__t33',t33)
                        p1.set('transformationMatrix__t34',t34)
			
			p1.setMicrograph(m1)
			emxData_p.addObject(p1)
		
			param = param_list.readline()
			j = j + 1			

	
	for _object in emxData_m:
		emxData.addObject(_object)

	for _object in emxData_p:
		emxData.addObject(_object)

	xmlMapper.writeEMXFile(emx_file)	

if __name__ == "__main__":
        main()
        sys.exit()
	
