Extraction of ECMWF data for the trajectory model FLEXTRA and the particle dispersion model FLEXPART

NEW UPDATE

If you need input data for FLEXTRA only, you can leave away step 1

Warning: In order to use these routines, it is necessary to have a user account on ecgate1.ecmwf.int

For more information: send email to Paul James


The tar file ecdocs.tar

All routines and documentation needed are enclosed in the file ecdocs.tar:

Step 1: Retrieval and de-accumulation of flux data

The programs can be found in the directory ./fluxjobs_gen and ./fluxjobs_era (for ERA-15 data) .

The following flux data needed by FLEXPART are retrieved on a user-selectable geographical grid (lamda-phi):

Abreviation Full name Unit Field code
LSP Large Scale Precipitation mm/hr 142
CP Convective Precipitation mm/hr 143
SSHF Surf. Sens. Heat Flux W/m2 146
EWSS East-West Surf. Stress N/m2 180
NSSS North-South Surf. Stress N/m2 181
SSR Surf. Solar Radiation W/m2 176

Flux data are not available in the ECMWF analyses, but they can be obtained (in accumulated form) from the ECMWF short-term forecasts (first guess) stored in three-hourly intervals (see ECMWF MARS database manual).

For de-accumulation, the following strategy is followed:

If you want to work with these routines, you should read the following instructions:

Step 2: Retrieval of other meteorological fields, calculation of the vertical velocity component d(eta)/dt

The programs can be found in the directories  ./job_era  ./job_pre99  ./job_y2k

This step has been fully automated. The user just has to fill in a CONTROL file, providing information like beginning date, ending date, time interval, number of levels to extract, scratch directory for (optional) flux data input, ecfile directory for output. Afterwards, the execution of a.out creates a batch script and sends it to batch queue ecgate1.normal.

The batch file is saved under /tmp/job.prepo.bak.uid. After execution, a job log file in ~uid/waitqueue named METEprepXXXXXX with XXXXXX being a system-generated, unique number, is created. Beginning and end of execution are anounced to the user via e-mail.

The batch-job works as follows:

To use these routines, please follow the instructions step by step: Technical comments:

Update August 1998 - Retrieval of additional parameters, 3 prepared domains

In August 1998, the retrieval routines (Step 2) were revised to get additional quantities potentially useful for particle modelling, air quality modelling and mesoscale meteorological modelling. The following additional data are now provided:
Abreviation Field code Name
SDOR 160 Std. Deviaton of Orography
VEG 199 Percentage of vegetation
SR 173 Surface roughness
LCC 186 Low cloud cover
MCC 187 Medium cloud cover
HCC 188 High cloud cover
SKT 235 Skin temperature
STL1 139 Soil temperature level 1
SWL1 140 Soil wetness level 1

Please notice: The new routines just work for ECMWF data from 04/04/1995 0 UTC onwards. If you want to extract older ECMWF data, use older routines written by Gerhard Wotawa.

If you want to use these old routines, just download them.
 

Handling of wind components at the poles

At the south Pole (phi=-90 DEG) and the north Pole (phi=90 DEG), all scalar quantities (e.g., temperature) are set to the same value for -178<= lambda <= 180 DEG (i=1,NX) on the lambda-phi grid. The wind components u, v are handled according to the following convention:

North Pole:

u(lambda)=-ff sin(lambda+dd)
v(lambda)=-ff cos(lambda+dd)

 

 

where lambda are the geographical longitude (-pi to + pi), ff the wind speed at the pole, dd the wind direction at the pole. Both, ff and dd, are constant for the Pole.

The wind direction is calculated as follows:

dd(lambda)=atan(u(lambda)/v(lambda))-lambda [if v(lambda)<0]
dd(lambda)=pi+atan(u(lambda)/v(lambda))-lambda [else]

 

 

For example: A wind direction of pi (180 DEG) means that the wind at the pole blows along the Greenwich meridian towards the pole. If the wind direction is 0, the wind blows the other way round.

South Pole:

u(lambda)=+ff sin(lambda-dd)
v(lambda)=-ff cos(lambda-dd)

 

 

The wind direction is calculated as follows:

dd(lambda)=atan(u(lambda)/v(lambda))+lambda [if v(lambda)<0]
dd(lambda)=pi+atan(u(lambda)/v(lambda))+lambda [else]

 

 

For a better understandig, I provide a little program where wind at the south and north pole is transformed to speed and direction and subsequently transformed back according to the pole convention used by ECMWF:

* SOUTH POLE 
 
      DO 666 I=1,NX
         LLON=-178+FLOAT(I-1)*DX
         LLON1=LLON*PI/180.
         WSPEED=SQRT((U3D(I,1,1)**2)+(V3D(I,1,1)**2))
         IF(V3D(I,1,1).lt.0.) THEN
           DDPOL=ATAN(U3D(I,1,1)/V3D(I,1,1))+LLON1
         ELSE
           DDPOL=PI+ATAN(U3D(I,1,1)/V3D(I,1,1))+LLON1
         ENDIF
         IF(DDPOL.LT.0.) DDPOL=2.0*pi+DDPOL
         IF(DDPOL.GT.2.0*PI) DDPOL=DDPOL-2.0*pi
         UNEW=+WSPEED*SIN(LLON1-DDPOL)
         VNEW=-WSPEED*COS(LLON1-DDPOL)
666      WRITE(999,778) I,LLON,U3D(I,1,1),V3D(I,1,1),WSPEED,
     &         UNEW,VNEW

*NORTH POLE 

            DO 777 I=1,NX
         LLON=-178+FLOAT(I-1)*DX
         LLON1=LLON*PI/180.
         WSPEED=SQRT((U3D(I,NY,1)**2)+(V3D(I,NY,1)**2))
         IF(V3D(I,NY,1).lt.0.) THEN
           DDPOL=ATAN(U3D(I,NY,1)/V3D(I,NY,1))-LLON1
         ELSE
           DDPOL=PI+ATAN(U3D(I,NY,1)/V3D(I,NY,1))-LLON1
         ENDIF
         IF(DDPOL.LT.0.) DDPOL=2.0*pi+DDPOL
         IF(DDPOL.GT.2.0*PI) DDPOL=DDPOL-2.0*pi
         UNEW=-WSPEED*SIN(LLON1+DDPOL)
         VNEW=-WSPEED*COS(LLON1+DDPOL)
777      WRITE(999,778) I,LLON,U3D(I,NY,1),V3D(I,NY,1),WSPEED,
     &         UNEW,VNEW
778      FORMAT(i3,1x,f5.0,1x,5f7.2)
Important remark: A spectral truncation of T106 is default in the ECMWF mars system for the retrieval of 2x2 DEG lambda-phi data out of the Sperical Harmonics (SH) data. However, as I have tested, this truncation is insufficient to obtain accurate wind data at the poles. If you use the latest (6.1.1999) version of my routines, this problem is fixed. If not, make sure that U,V wind components are retrieved at least with T213 or T319. This problem is serious. Wind data at the poles are definitely unusable with the lower truncation.

Bug reports:

(1) Date: 12.1.1999
Modification: File makeinclude.f in directory job.tra1995+ modified
Please note: Values of etapoint and log surface pressure are wrong in the global 2x2 grid as well as in the 0.4 degrees grid prior to this modification.

Restrictions

This software can be used "as is" by all people having user accounts at the ECMWF workstation ecgate1.ecmwf.int, without any fee or cost. Any warranty is disclaimed and no support is given. This software is part of the FLEXTRA and FLEXPART program package. The same regulations as accepted by using FLEXTRA and FLEXPART do apply. Please keep in mind the data usage regulations of the European Centre for Medium Range Weather Forecasts and of your national weather service. Any use of this program, which is against these regulations, is prohibited.

You are kindly invited to report problems and bugs to me.





Appendix: Plot of sample fields

Last update: 7 August 2000  by Paul James