fout = "cluster_model.dat"


N = 10000L


pos = fltarr(3, N)
vel = fltarr(3, N)
id = lindgen(N)+1

mp =  1.0  ;;;; particle mass



;;; let's fill in the coordinates -- just some random numbers here. Replace.

for i=0L, N-1 do begin

   pos(0,i) = randomu(seed)
   pos(1,i) = randomu(seed)
   pos(2,i) = randomu(seed)

endfor

;;; let's fill in the velocities  -- just some random numbers here. Replace.

for i=0L, N-1 do begin

   vel(0,i) = randomu(seed)
   vel(1,i) = randomu(seed)
   vel(2,i) = randomu(seed)

endfor





;;; Now write this as a GADGET2 file 

;;; first, create the variables in the header

npart=lonarr(6)
massarr=dblarr(6)
time=0.0D
redshift=0.0D
flag_sfr=0L
flag_feedback=0L
npartTotal=lonarr(6)
flag_cooling=0L
num_files = 1L
BoxSize=double(0)
Omega0=0.0D
OmegaLambda=0.0D
HubbleParam=0.0D
flag_stellarage=0L
flag_metals=0L
highwork=lonarr(6)
flag_entropy=0L
flag_double=0L
flag_lpt=0L
factor=0.0
la=intarr(48/2)

;;;; fill in header information

npart(1)      = N       ;;; we use the "halo" particle type 1
npartTotal(1) = N
massarr(1)    = mp      ;;; all particles have this mass


;;;;; write file

openw, 1, fout, /f77_unformatted

writeu,1,npart,massarr,time,redshift,flag_sfr,flag_feedback,npartTotal, $
       flag_cooling, num_files, BoxSize, Omega0, OmegaLambda, HubbleParam, flag_stellarage, flag_metals, highwork,$
       flag_entropy, flag_double, flag_lpt, factor, la

writeu,1,pos
writeu,1,vel
writeu,1,id

close,1


end