GRIB är ett binärt dataformat för data på ett beräkningsbart gridnät. Filen innehåller mestadels kod som kombineras med externa metadata för att rätt kunna läsa ut och tolka datasetet i GRIB-filen. Här finner du en översiktlig beskrivning av GRIB-formatet.
Nedladdningsbara meteorologiska modeller har ett roterat latitud/longitud gridnät medan oceanografiska modeller har ett reguljärt latitud/longitud gridnät.
För att kunna läsa data på GRIB-format behöver du:
- Program för att öppna och avkoda GRIB-fil
- Översätta möjliga roterade koordinater
- Rotera möjliga roterade vindar
- Kodtabell för att tolka parameternummer
Program för läsande av GRIB-fil
För att kunna läsa GRIB-packad data behövs program för att packa upp GRIB-filen och tolka data. Vid SMHI används mestadels program framtagna av ECMWF (http://www.ecmwf.int).
Dokumentation av den senaste varianten av programpaket (GRIB-API) samt möjlighet att ladda ned källkoden finns på ECMWF:s hemsida för GRIB-API,
https://software.ecmwf.int/wiki/display/GRIB/Home.
Andra program för läsning, skrivning och visualisering av data på GRIB-format finns att tillgå.
Roterade koordinater
De meteorologiska modeller som kan hämtas som ”Öppna data” har data lagrat på ett roterat gridnät. Roterade koordinater innebär att man flyttar sydpolen från lat/lon -90/0. Skälet till rotationen är att man genom att flytta ekvatorn till Sveriges latitud får mindre skillnad i avståndet mellan gridpunkter i norra o södra delen av modellområdet än vid oroterat gridnät.
I de fall då gridnätet är roterat kan koordinaterna i gridnätet behöva beräknas om till reguljär latitud och longitud. Detta kan göras med hjälp av Fortran-rutinen REGROT, se nedan, eller de ekvationer den innehåller. All indata som behövs för REGROT finns som metadata i GRIB-filen.
Rotera vinden
Vinden beskrivs i de nedladdningsbara GRIB-filerna som U/V (ost/nord) relaterat till det underliggande gridnätet. Är gridnätet roterat är även vinden roterad och för att bli korrekt i ett reguljärt gridnät behöver en omräkning ske.
Denna omräkning kan ske med hjälp av Fortran-rutinen TURNWI, se nedan, eller de ekvationer den innehåller.
Parameter-koder
Metadata för att beskriva vad en parameterkod betyder finns beskrivet i en extern fil. Olika tabeller översätter olika modellers parametrar. Vilken tabellversion som används i en fil anges av metadata-nyckeln table2Version i sektion 1 i GRIB-filen.
Kodtabeller för de modeller som finns nedladdningsbara kan laddas ner som en textfil; Kodlista Arome (3 kB, txt), Kodlista Arome (3 kB, txt), Kodlista Mesan (1 kB, txt), Kodlista PMP (2 kB, txt), Kodlista Hiromb (1 kB, txt).
Gridnät för vindberäkningar
Det finns en skillnad i hur beräkningar av vinden sker på olika nivåer. U och V beräknas på olika typer av ARAKAWA-grid beroende på typ av nivå.
Vinden på 10-meters nivå beräknas enligt ARAKAWA C medan vinden på modellnivåer (indicatorOfTypeOfLevel=109) beräknas enligt ARAKAWA A.
Exempel på rotation och område
En modells hörnkoordinater (beräknat i modellens projektion) och rotation finns beskrivet i form av metadata i GRIB-meddelandet. I GRIB edition 1 anges detta (värden tagna från PMP) bl.a. med metadata ur sektion 2:
latitudeOfFirstGridPoint = -6749 [milligrader]
longitudeOfFirstGridPoint = -7775 [milligrader]
latitudeOfSouthernPole = -30000 [milligrader]
longitudeOfSouthernPole = 15000 [milligrader]
De två första värdena anger i PMP:s fall sydvästra hörnet i modellens projektion (roterad latitud/longitud i PMP). De två sista talar om till vilka koordinater sydpolen har flyttats till i samband med rotationen.
Omräkning mellan roterat och reguljär projektion är komplicerat. Som en förenklad sett att förklara områdets storlek kan man tänka sig att skärningen mellan nollmeridianen och ekvatorn genom rotation flyttats till 15° ost och 60° nord. Utifrån denna punkt kan man grovt ange sydvästra hörnet av området till 6,749° syd och 7,775° väst.
Nordöstra hörnet kan sedan grovt beräknat med hjälp av avståndet mellan gridpunkterna och antalet punkter i x och y-led.
Kod för regrot
CALL REGROT(ORTLON,ORTLAT,
+ ALON_R,ALAT_R,
+ 1,1,1,1,ALONSP,ALATSP,1)
c ALON_R och ALAT_R är positionen i det roterade gridet
YY = 1. + (ALAT_R - ALATS)/DTHETA
XX = 1. + (ALON_R - ALONW)/DLAMDA
CALL INTER(NLON,NLAT,XX,YY,GRDVAL,PROG)
BIAS = GRDVAL-OBSVAL
WRITE(*,'(2f7.2,3f10.3)') ORTLAT,ORTLON,BIAS,
Fortrandkod för regrot
subroutine regrot(pxreg,pyreg,pxrot,pyrot,kxdim,kydim,kx,ky,
+ pxcen,pycen,kcall)
c
implicit none
c
c-------------------------------------------------------------------------------
c
c* conversion between regular and rotated spherical coordinates.
c*
C* pxreg longitudes of the regular coordinates
c* pyreg latitudes of the regular coordinates
c* pxrot longitudes of the rotated coordinates
c* pyrot latitudes of the rotated coordinates
c* all coordinates given in degrees N (negative values for S)
c* and degrees E (negative values for W)
c* kxdim dimension of the gridpoint fields in the x-direction
c* kydim dimension of the gridpoint fieLds in the y-direction
c* kx number of gridpoint in the x-direction
c* ky number of gridpoints in the y-direction
c* pxcen regular longitude of the south pole of the rotated grid
c* pycen regular latitude of the south pole of the rotated grid
c*
C* kcall=-1: find regular as functions of rotated coordinates.
c* kcall= 1: find rotated as functions of regular coordinates.
c*
c* j.e. haugen hirlam june -92
c
c-------------------------------------------------------------------------------
c
integer kxdim,kydim,kx,ky,kcall
real pxreg(kxdim,kydim),pyreg(kxdim,kydim),
+ pxrot(kxdim,kydim),pyrot(kxdim,kydim),
+ pxcen,pycen
c
c-------------------------------------------------------------------------------
c
real pi,zrad,zsycen,zcycen,zxmxc,zsxmxc,zcxmxc,zsyreg,zcyreg,
x zsyrot,zcyrot,zcxrot,zsxrot,zradi
integer jy,jx
c
c-------------------------------------------------------------------------------
c
pi = 4.*atan(1.)
zrad = pi/180.
zradi = 1./zrad
zsycen = sin(zrad*(pycen+90.))
zcycen = cos(zrad*(pycen+90.))
c
if (kcall.eq.1) then
c
do jy = 1,ky
do jx = 1,kx
c
zxmxc = zrad*(pxreg(jx,jy) - pxcen)
zsxmxc = sin(zxmxc)
zcxmxc = cos(zxmxc)
zsyreg = sin(zrad*pyreg(jx,jy))
zcyreg = cos(zrad*pyreg(jx,jy))
zsyrot = zcycen*zsyreg - zsycen*zcyreg*zcxmxc
zsyrot = amax1(zsyrot,-1.0)
zsyrot = amin1(zsyrot,+1.0)
c
pyrot(jx,jy) = asin(zsyrot)*zradi
c
zcyrot = cos(pyrot(jx,jy)*zrad)
zcxrot = (zcycen*zcyreg*zcxmxc +
+ zsycen*zsyreg)/zcyrot
zcxrot = amax1(zcxrot,-1.0)
zcxrot = amin1(zcxrot,+1.0)
zsxrot = zcyreg*zsxmxc/zcyrot
c
pxrot(jx,jy) = acos(zcxrot)*zradi
c
if (zsxrot.lt.0.0) pxrot(jx,jy) = -pxrot(jx,jy)
c
enddo
enddo
c
elseif (kcall.eq.-1) then
c
do jy = 1,ky
do jx = 1,kx
c
zsxrot = sin(zrad*pxrot(jx,jy))
zcxrot = cos(zrad*pxrot(jx,jy))
zsyrot = sin(zrad*pyrot(jx,jy))
zcyrot = cos(zrad*pyrot(jx,jy))
zsyreg = zcycen*zsyrot + zsycen*zcyrot*zcxrot
zsyreg = amax1(zsyreg,-1.0)
zsyreg = amin1(zsyreg,+1.0)
c
pyreg(jx,jy) = asin(zsyreg)*zradi
c
zcyreg = cos(pyreg(jx,jy)*zrad)
zcxmxc = (zcycen*zcyrot*zcxrot -
+ zsycen*zsyrot)/zcyreg
zcxmxc = amax1(zcxmxc,-1.0)
zcxmxc = amin1(zcxmxc,+1.0)
zsxmxc = zcyrot*zsxrot/zcyreg
zxmxc = acos(zcxmxc)*zradi
if (zsxmxc.lt.0.0) zxmxc = -zxmxc
c
pxreg(jx,jy) = zxmxc + pxcen
c
enddo
enddo
c
else
write(6,'(1x,''invalid kcall in regrot'')')
stop
endif
c
return
end
Fortrandkod för TURNWI
subroutine turnwi(puarg,pvarg,pures,pvres,pa,pb,pc,pd,
+ pxreg,pyreg,pxrot,pyrot,
+ pxcen,pycen,kcall)
c
c* turn horizontal velocity components between regular and
c* rotated spherical coordinates.
c* kcall=-2 or 2: preparation of coefficients.
c* kcall=-1 or 1: multiplication after preparations.
c*
c* j.e. haugen hirlam june -92
c
real puarg,pvarg,
+ pures,pvres,
+ pa,pb,
+ pc,pd,
+ pxreg,pyreg,
+ pxrot,pyrot
c
if (kcall.eq.1 .or. kcall.eq.-1) then
c
c* multiplication between regular and rotated spherical grid
c
pures = pa*puarg + pb*pvarg
pvres = pc*puarg + pd*pvarg
c
elseif (kcall.eq.2) then
c
c* precalculate matrixes from regular to rotated grid
c
zsyc = sin(pycen)
zcyc = cos(pycen)
c
zsxreg = sin(pxreg)
zcxreg = cos(pxreg)
zsyreg = sin(pyreg)
zcyreg = cos(pyreg)
c
zxmxc = pxreg - pxcen
zsxmxc = sin(zxmxc)
zcxmxc = cos(zxmxc)
c
zsxrot = sin(pxrot)
zcxrot = cos(pxrot)
zsyrot = sin(pyrot)
zcyrot = cos(pyrot)
c
pa = zcyc*zsxmxc*zsxrot + zcxmxc*zcxrot
pb = zcyc*zcxmxc*zsyreg*zsxrot - zsyc*zcyreg*zsxrot -
+ zsxmxc*zsyreg*zcxrot
pc = zsyc*zsxmxc/zcyrot
pd = (zsyc*zcxmxc*zsyreg + zcyc*zcyreg)/zcyrot
c
elseif (kcall.eq.-2) then
c
c* precalculate matrixes from rotated to regular grid
c
zsyc = sin(pycen)
zcyc = cos(pycen)
c
zsxreg = sin(pxreg)
zcxreg = cos(pxreg)
zsyreg = sin(pyreg)
zcyreg = cos(pyreg)
c
zxmxc = pxreg - pxcen
zsxmxc = sin(zxmxc)
zcxmxc = cos(zxmxc)
c
zsxrot = sin(pxrot)
zcxrot = cos(pxrot)
zsyrot = sin(pyrot)
zcyrot = cos(pyrot)
c
pa = zcxmxc*zcxrot + zcyc*zsxmxc*zsxrot
pb = zcyc*zsxmxc*zcxrot*zsyrot + zsyc*zsxmxc*zcyrot -
+ zcxmxc*zsxrot*zsyrot
pc =-zsyc*zsxrot/zcyreg
pd = (zcyc*zcyrot - zsyc*zcxrot*zsyrot)/zcyreg
c
else
write(6,'(1x,''invalid kcall in turnwi'')')
stop
endif
c
return
end