import struct
import os
import os.path
import shutil
import math
import subprocess
import pickle
from numpy import * 
import numpy as numpy
import pickle
from matplotlib import pyplot
import csv
from ctypes import *
import struct

class atributes(Structure):
   _fields_ = [("pos", POINTER(c_float)),
               ("vel",POINTER( c_float)),
               ("masses",c_float),
               ("id",c_int),
               ("u",c_float)]

class particles(Structure):
    ngas=10000
    _fields_ = [('gas',atributes*ngas)]#,
                #('dark_matter',atributes*n_dm),
               # ('disk',atributes*ndisk),
               # ('bulge',atributes*nbulge),
              #  ('stars',atributes*nstars)]



class Create_Header(Structure):

   SixDoubles=c_double*6
   SixInts=c_int*6    
   OneDouble=c_double*1
   _fields_ = [("npart",SixInts),
               ("massarr",SixDoubles),
               ("time", c_float ),
               ("redshift",  c_float),
               ("flag_sfr", c_int),
               ("flag_feedback", c_int),
               ("npartTotal",SixInts),
               ("flag_cooling", c_int),
               ("num_files",c_int ),
               ("BoxSize", c_float),
               ("Omega0", c_float),
               ("OmegaLambda", c_float),
               ("HubbleParam",c_float),
               ("flag_stellarage",c_int ),
               ("flag_metals",c_int ),
               ("highwork",SixInts),
               ("flag_entropy",c_int ),
               ("flag_double",c_int ),
               ("flag_lpt",c_int),
               ("factor", c_float),
               ("bytes",c_byte*(256- 6*4- 6*8- 2*4- 2*4- 6*4- 2*4 - 4*4- 2*4- 6*4- 3*4 -4))]



def Random_Distribution(ngas):
    theta = []
    r = []
    w=[]
    for i in range(ngas):
        u1 = random.random()
        u2 = random.random()
        theta.append(u1*2.0*numpy.pi)
        r.append(numpy.sqrt(-2.0*numpy.log(u2)))
    x=r*numpy.cos(theta)
    y=r*numpy.sin(theta)  
    return x,y

def plotting(x,y):
    fig = pyplot.figure()
    pyplot.minorticks_on()
    ax1 = fig.add_subplot(111)
    ax1.minorticks_on()
    ax1.plot(x,y, 'ro')
    ax1.set_xlabel('$x$ $(R_{Sun})$',fontsize=15, weight='bold')
    ax1.set_ylabel(r'$y$ $(R_{Sun})$',fontsize=15, weight='bold')
    pyplot.savefig("testing_one.png")


if __name__ == "__main__":


    RSun=6.9598E10 # cm
    MSun=1.9892E33 # g

    mass_star=1.0*MSun
    radius_star =1.0*RSun

    ngas=10000
    boxsize=2*radius_star
    x,y=Random_Distribution(ngas)
   # plotting(x,y)



#Create header..........................................

    SixDoubles=c_double*6
    SixInts=c_int*6

    header=Create_Header()
    header.npart= SixInts(ngas,0,0,0,0,0) 
    header.massarr=SixDoubles(0.0,0.0,0.0,0.0,0.0,0.0) 
    header.time=0.0
    header.redshift= 0.0 
    header.flag_sfr= 0 
    header.flag_feedback= 0 
    header.npartTotal= SixInts(ngas,0,0,0,0,0) 
    header.flag_cooling=0  
    header.num_files=1  
    header.BoxSize= boxsize 
    header.Omega0=0.0  
    header.OmegaLambda =0.0  
    header.HubbleParam =0.0  
    header.flag_stellarage=0 
    header.flag_metals=0 
    header.highwork=SixInts(0,0,0,0,0,0) 
    header.flag_entropy=0 
    header.flag_double=0 
    header.flag_lpt=0 
    header.factor=0.0 



    print 'how many?',sizeof(header)
    print 'writing header in file'

    f = open("header.dat","wb")
    p = struct.pack("iiiiiiffffffffiiiiiiiiiiffffiiiiiiiiiiif",header.npart[0],header.npart[1],header.npart[2],
                        header.npart[3],header.npart[4],header.npart[5],\
                        header.massarr[0],header.massarr[1],header.massarr[2],header.massarr[3],header.massarr[4],header.massarr[5],\
                        header.time,header.redshift,header.flag_sfr,header.flag_feedback,\
                        header.npartTotal[0],header.npartTotal[1],header.npartTotal[2],header.npartTotal[3],header.npartTotal[4],\
                        header.npartTotal[5], header.flag_cooling, header.num_files, header.BoxSize, header.Omega0,header.OmegaLambda,header.HubbleParam,\
                        header.flag_stellarage, header.flag_metals,header.highwork[0],header.highwork[1],header.highwork[2],header.highwork[3],\
                        header.highwork[4],header.highwork[5],header.flag_entropy, header.flag_double, header.flag_lpt, header.factor)
   
    f.write(p)

    star1=particles()
    atributes.pos=(c_float*3)(0.0,0.0,0.0)
    atributes.vel=(c_float*3)(0.0,0.0,0.0)
    atributes.masses=(c_float*1)(0.0)
    atributes.id=(c_int*1)(0)

    positions=range(3)
    velocities=range(3)

    for i in range(ngas):
        star1.gas[i].pos[0]= x[i]
        star1.gas[i].pos[1]= y[i]
        star1.gas[i].pos[2]= 0.0
        for k in range(3):
            positions[k]=star1.gas[i].pos[k]
        f.write(struct.pack('f'*(len(positions)), *positions))

        star1.gas[i].vel[0]= 0.0
        star1.gas[i].vel[1]= 0.0
        star1.gas[i].vel[2]= 0.0
        for k in range(3):
            velocities[k]=star1.gas[i].vel[k]
        f.write(struct.pack('f'*(len(velocities)), *velocities))

        star1.gas[i].id[0]=(i+1) 
        f.write(struct.pack('i',star1.gas[i].id[0] )) 

        star1.gas[i].masses[0]=(mass_star/ngas)
        f.write(struct.pack( 'f',star1.gas[i].masses[0]))

        star1.gas[i].u=(0.0)
        f.write(struct.pack( 'f',star1.gas[i].u))

    f.close()


    g=open('header.dat','rb')
    q=header
    print 'how many?',sizeof(header)
    print header.time