Index: README_PARALLEL.md
===================================================================
--- README_PARALLEL.md	(revision 1d207bb19e353f44d1fa73e9527e91a88652a8bd)
+++ README_PARALLEL.md	(revision 5f42c27a3ad632d81e1e6ecd2f81d821b0c11373)
@@ -1,4 +1,4 @@
 
-		FLEXPART VERSION 9.2.0 (MPI)
+		FLEXPART VERSION 10.0 beta (MPI)
 
 Description
@@ -62,5 +62,5 @@
   'numpart' are complemented by variables 'maxpart_mpi' and 'numpart_mpi'
   which are the run-time determined number of particles per process, i.e,
-  maxpart_mpi = maxpart/[number of processes]. The variable 'numpart' 
+  maxpart_mpi = maxpart/np, where np are the number of processes. The variable 'numpart' 
   is still used in the code, but redefined to mean 'number of particles
   per MPI process'
@@ -79,5 +79,5 @@
   But it is also possible that the
   program will run even faster if the 4th process is participating in 
-  the calculation of particle trajectories. This will largely depend on
+  the calculation of particle trajectories instead. This will largely depend on
   the problem size (total number of particles in the simulation, resolution
   of grids etc) and hardware being used (disk speed/buffering, memory
@@ -120,5 +120,5 @@
   FLEXPART-MPI is run with no (or just a few) particles, no performance 
   improvement is possible. In this case, most processing time is spent
-  in the 'getfields'-routine (ECMWF).
+  in the 'getfields'-routine.
 
   A) Running without dedicated reader process
@@ -173,5 +173,5 @@
     * Nested wind fields
 
- -The following will probably/possibly not work (untested/under developement): 
+ -The following will most probably not work (untested/under developement): 
 
     * Backward runs
Index: options/COMMAND
===================================================================
--- options/COMMAND	(revision 5f9d14aa526df192efa539f87fc409fdeeeee06c)
+++ options/COMMAND	(revision 5f42c27a3ad632d81e1e6ecd2f81d821b0c11373)
@@ -12,5 +12,5 @@
  CTL= -5.0000000    ,
  IFINE=          4,
- IOUT=          1,
+ IOUT=          9,
  IPOUT=          0,
  LSUBGRID=          1,
@@ -28,3 +28,4 @@
  SURF_ONLY=          0,
  CBLFLAG=          0,
+ OHFIELDS_PATH=          "../../flexin/",
  /
Index: options/OUTGRID
===================================================================
--- options/OUTGRID	(revision 60403cdfd1dfc5def4c83aff4c9c55a236804c93)
+++ options/OUTGRID	(revision 5f42c27a3ad632d81e1e6ecd2f81d821b0c11373)
@@ -6,71 +6,151 @@
 ********************************************************************************
 
+1.  ------.----       4X,F11.4
+      -179.0000       GEOGRAFICAL LONGITUDE OF LOWER LEFT CORNER OF OUTPUT GRID
+    OUTLONLEFT        (left boundary of the first grid cell - not its centre)
 
-      -179.0000
+2.  ------.----       4X,F11.4
+       -90.0000       GEOGRAFICAL LATITUDE OF LOWER LEFT CORNER OF OUTPUT GRID
+    OUTLATLOWER       (lower boundary of the first grid cell - not its centre)
+
+3.  -----             4X,I5
+     720               NUMBER OF GRID POINTS IN X DIRECTION (= No. of cells + 1)
+    NUMXGRID
+
+4.  -----             4X,I5
+     360             NUMBER OF GRID POINTS IN Y DIRECTION (= No. of cells + 1)
+    NUMYGRID
+
+5.  ------.---        4X,F10.3
+         0.50        GRID DISTANCE IN X DIRECTION
+    DXOUTLON
+
+6.  ------.---        4X,F10.3
+         0.50        GRID DISTANCE IN Y DIRECTION
+    DYOUTLAT
+
+10. -----.-           4X, F7.1
+     1000.0
+    LEVEL 1           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     2000.0
+    LEVEL 2           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     3000.0
+    LEVEL 3           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     4000.0
+    LEVEL 4           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     5000.0
+    LEVEL 5           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     6000.0
+    LEVEL 6           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     7000.0
+    LEVEL 7           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     8000.0
+    LEVEL 8           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+     9000.0
+    LEVEL 9           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    10000.0
+    LEVEL10           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    11000.0
+    LEVEL11           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    12000.0
+    LEVEL12           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    13000.0
+    LEVEL13           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    14000.0
+    LEVEL14           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    15000.0
+    LEVEL15           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    16000.0
+    LEVEL16           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    17000.0
+    LEVEL17           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    18000.0
+    LEVEL18           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    19000.0
+    LEVEL19           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    20000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    21000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    22000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    23000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    24000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    25000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    26000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    27000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    28000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    29000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    30000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
+
+10. -----.-           4X, F7.1
+    40000.0
+    LEVEL20           HEIGHT OF LEVEL (UPPER BOUNDARY)
 
 
-
-       -90.0000
-
-
-
-      360
-
-
-
-      180
-
-
-
-         1.000
-
-
-
-         1.000
-
-
-
-      400.0
-
-
-
-     1000.0
-
-
-
-     2000.0
-
-
-
-     3000.0
-
-
-
-     5000.0
-
-
-
-     7000.0
-
-
-
-     9000.0
-
-
-
-    12000.0
-
-
-
-    15000.0
-
-
-
-    20000.0
-
-
-
-    50000.0
-
-
-
Index: src/makefile
===================================================================
--- src/makefile	(revision 5f9d14aa526df192efa539f87fc409fdeeeee06c)
+++ src/makefile	(revision 5f42c27a3ad632d81e1e6ecd2f81d821b0c11373)
@@ -34,4 +34,5 @@
 FLEXPART-GFS-MPI 	= FP_gfs_MPI
 
+
 ## COMPILERS
 F90       = /usr/bin/gfortran
@@ -50,17 +51,18 @@
 
 ## 1) System libraries at NILU
-# INCPATH1 = /usr/include
-INCPATH1 = /xnilu_wrk/flex_wrk/bin64/grib_api/include
+INCPATH1 = /usr/include
+#INCPATH1 = /xnilu_wrk/flex_wrk/bin64/grib_api/include
 LIBPATH1 = /xnilu_wrk/flex_wrk/bin64/grib_api/lib
-# LIBPATH2 = /usr/lib/x86_64-linux-gnu
+#LIBPATH2 = /usr/lib/x86_64-linux-gnu
 
 ## 2) Home-made libraries
-#INCPATH2  = ${HOME}/include
-INCPATH2  = /homevip/flexpart/include
-LIBPATH2 = /homevip/flexpart/lib
-#INCPATH2 = ${HOME}/lib	
-
-FFLAGS   = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV) -g -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV) -mtune=native $(FUSER) # -march=native
-DBGFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV_DBG) -g3 -ggdb3 -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV_DBG) -fbacktrace -Warray-bounds  -Wall -fcheck=all $(FUSER)  # -ffpe-trap=invalid,overflow,denormal,underflow,zero -fdump-core
+# INCPATH1  = ${HOME}/include
+INCPATH2  = /homevip/flexpart/include/
+LIBPATH2 = /homevip/flexpart/lib/
+# LIBPATH1 = ${HOME}/lib
+
+FFLAGS   = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV) -g -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV) -mtune=native -fuse-linker-plugin $(FUSER) # -march=native
+
+DBGFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV_DBG) -g3 -ggdb3 -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV_DBG) -fbacktrace -Warray-bounds  -Wall -fcheck=all $(FUSER)  #  -ffpe-trap=invalid,overflow,denormal,underflow,zero -fdump-core
 
 LDFLAGS  = $(FFLAGS) -L$(LIBPATH1) -L$(LIBPATH2) $(LIBS)
@@ -188,4 +190,7 @@
 %.o: %.mod
 
+ecmwf: $(FLEXPART-ECMWF)
+ecmwf: FC := $(F90)
+
 ecmwf-mpi: $(FLEXPART-ECMWF-MPI)
 ecmwf-mpi: FC := $(MPIF90)
@@ -196,7 +201,4 @@
 ecmwf-mpi-dbg: FC := $(MPIF90)
 
-ecmwf: $(FLEXPART-ECMWF)
-ecmwf: FC := $(F90)
-
 gfs: $(FLEXPART-GFS)
 gfs: FC := $(F90)
@@ -207,4 +209,7 @@
 #all: $(FLEXPART-ECMWF-MPI)
 
+$(FLEXPART-ECMWF): $(MODOBJS) $(OBJECTS) $(OBJECTS_SERIAL) $(OBJECTS_ECMWF)
+	+$(FC) -o $@ $(MODOBJS) $(OBJECTS) $(OBJECTS_SERIAL) $(OBJECTS_ECMWF) $(LDFLAGS)
+
 $(FLEXPART-ECMWF-MPI): $(MODOBJS) $(MPI_MODOBJS) $(OBJECTS) $(OBJECTS_MPI) $(OBJECTS_ECMWF_MPI)
 	+$(FC) -o $@ $(MODOBJS) $(MPI_MODOBJS) $(OBJECTS) $(OBJECTS_MPI) \
@@ -216,7 +221,4 @@
 	+$(FC) -o $@ $(MODOBJS) $(MPI_MODOBJS) $(OBJECTS) $(OBJECTS_MPI) \
 	$(OBJECTS_ECMWF_MPI) $(LDFLAGS)
-
-$(FLEXPART-ECMWF): $(MODOBJS) $(OBJECTS) $(OBJECTS_SERIAL) $(OBJECTS_ECMWF)
-	+$(FC) -o $@ $(MODOBJS) $(OBJECTS) $(OBJECTS_SERIAL) $(OBJECTS_ECMWF) $(LDFLAGS)
 
 $(FLEXPART-GFS): $(MODOBJS) $(OBJECTS) $(OBJECTS_SERIAL) $(OBJECTS_GFS)
Index: c/random.f90
===================================================================
--- src/random.f90	(revision e200b7adb07d4dfa7bf42f7cb43e0c85b4dbf894)
+++ 	(revision )
@@ -1,154 +1,0 @@
-!**********************************************************************
-! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
-! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
-! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
-!                                                                     *
-! This file is part of FLEXPART.                                      *
-!                                                                     *
-! FLEXPART is free software: you can redistribute it and/or modify    *
-! it under the terms of the GNU General Public License as published by*
-! the Free Software Foundation, either version 3 of the License, or   *
-! (at your option) any later version.                                 *
-!                                                                     *
-! FLEXPART is distributed in the hope that it will be useful,         *
-! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
-! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
-! GNU General Public License for more details.                        *
-!                                                                     *
-! You should have received a copy of the GNU General Public License   *
-! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
-!**********************************************************************
-
-!  Taken from Press et al., Numerical Recipes
-
-function ran1(idum)
-
-  implicit none
-
-  integer :: idum
-  real    :: ran1
-  integer,parameter :: ia=16807, im=2147483647, iq=127773, ir=2836
-  integer,parameter :: ntab=32, ndiv=1+(im-1)/ntab
-  real,parameter    :: am=1./im, eps=1.2e-7, rnmx=1.-eps
-  integer :: j, k
-  integer :: iv(ntab) = (/ (0,j=1,ntab) /)
-  integer :: iy=0
-
-  if (idum.le.0.or.iy.eq.0) then
-    idum=max(-idum,1)
-    do j=ntab+8,1,-1
-      k=idum/iq
-      idum=ia*(idum-k*iq)-ir*k
-      if (idum.lt.0) idum=idum+im
-      if (j.le.ntab) iv(j)=idum
-    enddo
-    iy=iv(1)
-  endif
-  k=idum/iq
-  idum=ia*(idum-k*iq)-ir*k
-  if (idum.lt.0) idum=idum+im
-  j=1+iy/ndiv
-  iy=iv(j)
-  iv(j)=idum
-  ran1=min(am*iy,rnmx)
-end function ran1
-
-
-function gasdev(idum)
-
-  implicit none
-
-  integer :: idum
-  real    :: gasdev, fac, r, v1, v2
-  integer :: iset = 0
-  real    :: gset = 0.
-  real, external :: ran3
-
-  if (iset.eq.0) then
-1   v1=2.*ran3(idum)-1.
-    v2=2.*ran3(idum)-1.
-    r=v1**2+v2**2
-    if(r.ge.1.0 .or. r.eq.0.0) go to 1
-    fac=sqrt(-2.*log(r)/r)
-    gset=v1*fac
-    gasdev=v2*fac
-    iset=1
-  else
-    gasdev=gset
-    iset=0
-  endif
-end function gasdev
-
-
-subroutine gasdev1(idum,random1,random2)
-
-  implicit none
-
-  integer :: idum
-  real :: random1, random2, fac, v1, v2, r
-  real, external :: ran3
-
-1   v1=2.*ran3(idum)-1.
-  v2=2.*ran3(idum)-1.
-  r=v1**2+v2**2
-  if(r.ge.1.0 .or. r.eq.0.0) go to 1
-  fac=sqrt(-2.*log(r)/r)
-  random1=v1*fac
-  random2=v2*fac
-  ! Limit the random numbers to lie within the interval -3 and +3
-  !**************************************************************
-  if (random1.lt.-3.) random1=-3.
-  if (random2.lt.-3.) random2=-3.
-  if (random1.gt.3.) random1=3.
-  if (random2.gt.3.) random2=3.
-end subroutine gasdev1
-
-
-function ran3(idum)
-
-  implicit none
-
-  integer :: idum
-  real :: ran3
-
-  integer,parameter :: mbig=1000000000, mseed=161803398, mz=0
-  real,parameter    :: fac=1./mbig
-  integer :: i,ii,inext,inextp,k
-  integer :: mj,mk,ma(55)
-
-  save inext,inextp,ma
-  integer :: iff = 0
-
-  if(idum.lt.0.or.iff.eq.0)then
-    iff=1
-    mj=mseed-iabs(idum)
-    mj=mod(mj,mbig)
-    ma(55)=mj
-    mk=1
-    do i=1,54
-      ii=mod(21*i,55)
-      ma(ii)=mk
-      mk=mj-mk
-      if(mk.lt.mz)mk=mk+mbig
-      mj=ma(ii)
-    end do
-    do k=1,4
-      do i=1,55
-        ma(i)=ma(i)-ma(1+mod(i+30,55))
-        if(ma(i).lt.mz)ma(i)=ma(i)+mbig
-      end do
-    end do
-    inext=0
-    inextp=31
-    idum=1
-  endif
-  inext=inext+1
-  if(inext.eq.56)inext=1
-  inextp=inextp+1
-  if(inextp.eq.56)inextp=1
-  mj=ma(inext)-ma(inextp)
-  if(mj.lt.mz)mj=mj+mbig
-  ma(inext)=mj
-  ran3=mj*fac
-end function ran3
-!  (C) Copr. 1986-92 Numerical Recipes Software US.
