#!/usr/bin/env python

import sys, commands, os
import random
import numpy as np

InitialPOS = open('DynamicsINPUTS/inicond_H2_FWHM/Output/ICON.txt', 'r').readlines()
Cell = [[6.9371841440755109,    0.0000000000000000,    0.0000000000000000],
[0.0000000000000000,    5.6641871349036883,    0.0000000000000000],
[0.0000000000000000,    0.0000000000000000,   23.7731707681984119]]
Nsurf = 24

for i in range(1, len(InitialPOS) + 1):
	if not os.path.exists('{0:04d}'.format( i ) ):
		os.makedirs('{0:04d}'.format( i ) )
	POSCAR = open('{0:04d}/POSCAR'.format( i ), 'w')

	POSCAR.write('H2 + Pt(211) - HD-NNP - PBEa57-DF2 - 1x2-4L 0 K - Ideal surface\n')
	POSCAR.write("   1.00000000000000\n")
	POSCAR.write("    {0:18.16f}   {1:18.16f}   {2:18.16f}\n".format( Cell[0][0], Cell[0][1], Cell[0][2] ))
	POSCAR.write("    {0:18.16f}   {1:18.16f}   {2:18.16f}\n".format( Cell[1][0], Cell[1][1], Cell[1][2] ))
	POSCAR.write("    {0:18.16f}   {1:18.16f}   {2:18.16f}\n".format( Cell[2][0], Cell[2][1], Cell[2][2] ))
	POSCAR.write("Pt H\n")
	POSCAR.write("{0:d} 2\n".format( Nsurf ))
	POSCAR.write("Selective dynamics\n")
	POSCAR.write("Cartesian\n")
        POSCAR.write("  0.0000000000000000  0.0000000000000000  0.0000000000000000  F  F  F\n")
        POSCAR.write("  0.0000000000000000  2.8320935674518437  0.0000000000000000  F  F  F\n")
        POSCAR.write("  4.6750941632945668  1.4160467837259203 -0.5701893268823994  F  F  F\n")
        POSCAR.write("  4.6750941632945668  4.2481403511777645 -0.5701893268823994  F  F  F\n")
        POSCAR.write("  2.4581638349332278  0.0000000000000000 -1.2474901842223600  F  F  F\n")
        POSCAR.write("  2.4581638349332278  2.8320935674518437 -1.2474901842223600  F  F  F\n")
        POSCAR.write("  0.2330891442173026  1.4160467837259203 -2.3253291402336389  F  F  F\n")
        POSCAR.write("  0.2330891442173026  4.2481403511777645 -2.3253291402336389  F  F  F\n")
        POSCAR.write("  4.8017532306534072  0.0000000000000000 -3.0344251394874706  F  F  F\n")
        POSCAR.write("  4.8017532306534072  2.8320935674518437 -3.0344251394874706  F  F  F\n")
        POSCAR.write("  2.5585435147651263  1.4160467837259203 -3.8715834947720698  F  F  F\n")
        POSCAR.write("  2.5585435147651263  4.2481403511777645 -3.8715834947720698  F  F  F\n")
        POSCAR.write("  0.2419552532915434  0.0000000000000000 -4.6840694026286398  F  F  F\n")
        POSCAR.write("  0.2419552532915434  2.8320935674518437 -4.6840694026286398  F  F  F\n")
        POSCAR.write("  4.8632846744275779  1.4160467837259203 -5.4881351966950698  F  F  F\n")
        POSCAR.write("  4.8632846744275779  4.2481403511777645 -5.4881351966950698  F  F  F\n")
        POSCAR.write("  2.5258135037382967  0.0000000000000000 -6.3439151094662609  F  F  F\n")
        POSCAR.write("  2.5258135037382967  2.8320935674518437 -6.3439151094662609  F  F  F\n")
        POSCAR.write("  0.2412346677687078  1.4160467837259203 -7.1380607846599107  F  F  F\n")
        POSCAR.write("  0.2412346677687078  4.2481403511777645 -7.1380607846599107  F  F  F\n")
        POSCAR.write("  4.8660240971523683  0.0000000000000000 -7.9556157764292497  F  F  F\n")
        POSCAR.write("  4.8660240971523683  2.8320935674518437 -7.9556157764292497  F  F  F\n")
        POSCAR.write("  2.5536293824605574  1.4160467837259203 -8.7731707681984101  F  F  F\n")
        POSCAR.write("  2.5536293824605574  4.2481403511777645 -8.7731707681984101  F  F  F\n")

	X_d = float( InitialPOS[i-1].split()[0] )
	Y_d = float( InitialPOS[i-1].split()[1] )
	H1 = np.dot( [X_d, Y_d, 0.], Cell )
	r = float( InitialPOS[i-1].split()[2] )
	theta = float( InitialPOS[i-1].split()[3] )
	phi = float( InitialPOS[i-1].split()[4] )
	H1[2] = -0.5*r*np.cos( np.radians( theta ) ) + 7.5
	H2 = [ r*np.sin( np.radians( theta ) )*np.cos( np.radians( phi ) ) + H1[0],
		r*np.sin( np.radians( theta ) )*np.sin( np.radians( phi ) ) + H1[1],
		0.5*r*np.cos( np.radians( theta ) ) + 7.5 ]

	POSCAR.write('{0:20.16f}{1:20.16f}{2:20.16f}  T  T  T\n'.format( H1[0], H1[1], H1[2] ) )
	POSCAR.write('{0:20.16f}{1:20.16f}{2:20.16f}  T  T  T\n'.format( H2[0], H2[1], H2[2] ) )

	POSCAR.write('\n')

	for j in range(Nsurf):
		POSCAR.write('    0.0000000000000000    0.0000000000000000    0.0000000000000000\n')

	POSCAR.write('{0:20.16f}{1:20.16f}{2:20.16f}\n'.format( float( InitialPOS[i-1].split()[5] ) / 1000., float( InitialPOS[i-1].split()[6] ) / 1000., float( InitialPOS[i-1].split()[7] ) / 1000. ) )
	POSCAR.write('{0:20.16f}{1:20.16f}{2:20.16f}\n'.format( float( InitialPOS[i-1].split()[8] ) / 1000., float( InitialPOS[i-1].split()[9] ) / 1000., float( InitialPOS[i-1].split()[10] ) / 1000. ) )
