Subversion Repositories lagranto.ecmwf

Rev

Rev 36 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 36 Rev 39
Line 15... Line 15...
15 
c     Set parameters
 15 
c     Set parameters
16 
c     --------------------------------------------------------------
 16 
c     --------------------------------------------------------------
17 
 
 17 
 
18 
c     Maximum number of starting positions
 18 
c     Maximum number of starting positions
19 
      integer           nmax
 19 
      integer           nmax
20 
      parameter         (nmax=4000000)
 20 
      parameter         (nmax=10000000)
21 
 
 21 
 
22 
c     Maximum number of model levels
 22 
c     Maximum number of model levels
23 
      integer           nlevmax
 23 
      integer           nlevmax
24 
      parameter         (nlevmax=200)
 24 
      parameter         (nlevmax=200)
25
25
Line 79... Line 79...
79 
      integer           start_n                        ! Number of coordinates
 79 
      integer           start_n                        ! Number of coordinates
80 
      real              start_lat(nmax)                ! Latitudes
 80 
      real              start_lat(nmax)                ! Latitudes
81 
      real              start_lon(nmax)                ! Longitudes
 81 
      real              start_lon(nmax)                ! Longitudes
82 
      real              start_lev(nmax)                ! Levels (depending on vertical unit)
 82 
      real              start_lev(nmax)                ! Levels (depending on vertical unit)
83 
      real              start_pre(nmax)                ! Level in hPa
 83 
      real              start_pre(nmax)                ! Level in hPa
- 
 
 84 
      integer           start_valid(nmax)              ! Flag whether starting height is valid
84 
      integer           reftime(6)                     ! Reference time
 85 
      integer           reftime(6)                     ! Reference time
85 
      character*80      vars(10)                       ! Name of output fields (time,lon,lat,p)
 86 
      character*80      vars(10)                       ! Name of output fields (time,lon,lat,p)
86 
      real,allocatable, dimension (:,:,:) :: tra       ! Trajectories (ntra,ntim,ncol)
 87 
      real,allocatable, dimension (:,:,:) :: tra       ! Trajectories (ntra,ntim,ncol)
87 
      real              latmin,latmax
 88 
      real              latmin,latmax
88 
      real              lonmin,lonmax
 89 
      real              lonmin,lonmax
Line 132... Line 133...
132 
      real,allocatable, dimension (:,:,:) :: fld0
 133 
      real,allocatable, dimension (:,:,:) :: fld0
133 
      real,allocatable, dimension (:,:,:) :: fld1
 134 
      real,allocatable, dimension (:,:,:) :: fld1
134 
      real,allocatable, dimension (:,:  ) :: sfc0
 135 
      real,allocatable, dimension (:,:  ) :: sfc0
135 
      real,allocatable, dimension (:,:)   :: sfc1
 136 
      real,allocatable, dimension (:,:)   :: sfc1
136 
      real,allocatable, dimension (:,:)   :: mask
 137 
      real,allocatable, dimension (:,:)   :: mask
- 
 
 138 
      integer           n_outside
137 
 
 139 
 
138 
c     Externals
140 
c     Externals
139 
      real              int_index3     ! 3d interpolation
 141 
      real              int_index3     ! 3d interpolation
140 
      external          int_index3
 142 
      external          int_index3
141 
      real              sdis           ! Speherical distance
 143 
      real              sdis           ! Speherical distance
Line 1140... Line 1142...
1140 
      start_n = 0
 1142 
      start_n = 0
1141 
      do i=1,vn
 1143 
      do i=1,vn
1142 
         do j=1,hn
 1144 
         do j=1,hn
1143 
 
 1145 
 
1144 
            start_n = start_n + 1
 1146 
            start_n = start_n + 1
- 
 
 1147 
            if ( start_n.ge.nmax ) then
- 
 
 1148 
               print*,' ERROR: maximum number of start points exceeded'
- 
 
 1149 
               stop
- 
 
 1150 
            endif
1145 
            start_lon(start_n) = lonlist(j)
 1151 
            start_lon(start_n) = lonlist(j)
1146 
            start_lat(start_n) = latlist(j)
 1152 
            start_lat(start_n) = latlist(j)
1147 
            start_lev(start_n) = levlist(i)
 1153 
            start_lev(start_n) = levlist(i)
1148 
 
 1154 
 
1149 
         enddo
 1155 
         enddo
Line 1300... Line 1306...
1300 
            start_pre(i) = mdv
 1306 
            start_pre(i) = mdv
1301 
         endif
 1307 
         endif
1302 
 
 1308 
 
1303 
      enddo
 1309 
      enddo
1304 
 
 1310 
 
1305 
c     Remove all starting points outside the domain
 1311 
c     Mark all starting points outside the domain
1306 
      i        = 1
 1312 
      n_outside = 0
1307 
      do while ( i.le.start_n )
 1313 
      do i=1,start_n
1308
1314
- 
 
 1315 
         start_valid(i) = 1
1309 
         if ( abs(start_pre(i)-mdv).lt.eps ) then
 1316 
         if ( abs(start_pre(i)-mdv).lt.eps ) then
1310
1317
- 
 
 1318 
            start_valid(i) = 0
1311 
            if ( vmode.ne.'grid') then
 1319 
            if ( vmode.ne.'grid' ) then
- 
 
 1320 
               if ( n_outside.lt.10 ) then
- 
 
 1321 
                  print*,'  Outside ',
1312 
             print*,'  Outside ', start_lon(i),start_lat(i),start_lev(i)
 1322 
     >                 start_lon(i),start_lat(i),start_lev(i)
- 
 
 1323 
                  n_outside = n_outside + 1
- 
 
 1324 
               elseif ( n_outside.eq.10 ) then
- 
 
 1325 
                  print*,'  Outside: more than 10 starting points'
- 
 
 1326 
                  n_outside = n_outside + 1
- 
 
 1327 
               else
- 
 
 1328 
                  n_outside = n_outside + 1
- 
 
 1329 
               endif
1313 
            endif
 1330 
            endif
1314
1331
1315 
            do j=i,start_n
 1332 
         endif
- 
 
 1333
- 
 
 1334 
      enddo
- 
 
 1335
1316 
               start_lon(j) = start_lon(j+1)
 1336
1317 
               start_lat(j) = start_lat(j+1)
 1337 
c     Only keep valid starting points
- 
 
 1338 
      j = 0
1318 
               start_pre(j) = start_pre(j+1)
 1339 
      do i=1,start_n
1319 
               start_lev(j) = start_lev(j+1)
 1340 
         if ( start_valid(i).eq.1 ) then
1320 
            enddo
 1341 
             j = j + 1
1321 
            start_n = start_n - 1
 1342 
             start_lon(j) = start_lon(i)
1322 
 
 - 
 
1323 
         else
 1343 
             start_lat(j) = start_lat(i)
1324
1344 
             start_pre(j) = start_pre(i)
1325 
            i = i + 1
 1345 
             start_lev(j) = start_lev(i)
1326 
 
 - 
 
1327 
         endif
 1346 
          endif
1328
- 
 
1329 
      enddo
 1347 
      enddo
- 
 
 1348 
      start_n = j
1330 
 
 1349 
 
1331 
c     Write some status information
 1350 
c     Write some status information
1332 
      latmin = start_lat(1)
 1351 
      latmin = start_lat(1)
1333 
      latmax = start_lat(1)
 1352 
      latmax = start_lat(1)
1334 
      lonmin = start_lon(1)
 1353 
      lonmin = start_lon(1)
Line 1346... Line 1365...
1346 
      print*,'  min(lat),max(lat) : ', latmin,latmax
 1365 
      print*,'  min(lat),max(lat) : ', latmin,latmax
1347 
      print*,'  min(lon),max(lon) : ', lonmin,lonmax
 1366 
      print*,'  min(lon),max(lon) : ', lonmin,lonmax
1348 
      print*,'  min(pre),max(pre) : ', premin,premax
 1367 
      print*,'  min(pre),max(pre) : ', premin,premax
1349 
      print*
 1368 
      print*
1350 
      print*,'  # starting points : ', start_n
 1369 
      print*,'  # starting points : ', start_n
- 
 
 1370 
      print*,'  # Ouside domain   : ', n_outside
1351 
      print*
 1371 
      print*
1352 
 
 1372 
 
1353
1373
1354 
c     ------------------------------------------------------------------
 1374 
c     ------------------------------------------------------------------
1355 
c     Write starting positions to output file
 1375 
c     Write starting positions to output file
Line 1363... Line 1383...
1363 
         vars(1)  ='time'
 1383 
         vars(1)  ='time'
1364 
         vars(2)  ='lon'
 1384 
         vars(2)  ='lon'
1365 
         vars(3)  ='lat'
 1385 
         vars(3)  ='lat'
1366 
         vars(4)  ='p'
 1386 
         vars(4)  ='p'
1367 
         vars(5)  ='level'
 1387 
         vars(5)  ='level'
- 
 
 1388 
 
1368 
         call wopen_tra(fid,ofile,start_n,1,5,reftime,vars,oformat)
 1389 
         call wopen_tra(fid,ofile,start_n,1,5,reftime,vars,oformat)
1369 
 
 1390 
 
1370 
         do i=1,start_n
 1391 
         do i=1,start_n
1371 
            tra(i,1,1) = 0.
 1392 
            tra(i,1,1) = 0.
1372 
            tra(i,1,2) = start_lon(i)
 1393 
            tra(i,1,2) = start_lon(i)