Subversion Repositories lagranto.wrf

Rev

Rev 15 | Rev 21 | Go to most recent revision | Show entire file | Regard whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 15 Rev 20
Line 21... Line 21...
21
c     Fixed parameters
21
c     Fixed parameters
22
      character*80 mapfile                          ! netCDF file with map transformation
22
      character*80 mapfile                          ! netCDF file with map transformation
23
      parameter    (mapfile ='wrfmap.nc')
23
      parameter    (mapfile ='wrfmap.nc')
24
      real         mdv                              ! missing data value
24
      real         mdv                              ! missing data value
25
      parameter    (mdv=-999.)
25
      parameter    (mdv=-999.)
-
 
26
      real         eps                              ! Numerical epsilon
-
 
27
      parameter    (eps= 0.001)
26
 
28
 
27
c     Numerical grid
29
c     Numerical grid
28
      real         xmin,xmax,ymin,ymax              ! Domain size
30
      real         xmin,xmax,ymin,ymax              ! Domain size
29
      real         dx,dy                            ! Horizontal resolution
31
      real         dx,dy                            ! Horizontal resolution
30
      integer      nx,ny,nz                         ! Grid dimensions
32
      integer      nx,ny,nz                         ! Grid dimensions
Line 517... Line 519...
517
      endif 
519
      endif 
518
 
520
 
519
c     Transform coordinates
521
c     Transform coordinates
520
      do i=1,npoints
522
      do i=1,npoints
521
 
523
 
-
 
524
          if ( ( abs(xx0(i)-mdv).gt.eps ).and.
-
 
525
     >         ( abs(yy0(i)-mdv).gt.eps )  )
-
 
526
     >    then
522
          rid=(xx0(i)-lonmin)/dlon+1.
527
             rid=(xx0(i)-lonmin)/dlon+1.
523
          rjd=(yy0(i)-latmin)/dlat+1.
528
             rjd=(yy0(i)-latmin)/dlat+1.
524
 
-
 
525
          xx0(i) = int_index3(x,nlon,nlat,1,rid,rjd,1.,mdv)
529
             xx0(i) = int_index3(x,nlon,nlat,1,rid,rjd,1.,mdv)
526
          yy0(i) = int_index3(y,nlon,nlat,1,rid,rjd,1.,mdv)
530
             yy0(i) = int_index3(y,nlon,nlat,1,rid,rjd,1.,mdv)
-
 
531
          else
-
 
532
             xx0(i) = mdv
-
 
533
             yy0(i) = mdv
-
 
534
          endif
527
 
535
 
528
      enddo
536
      enddo
529
 
537
 
530
c     Write output to list
538
c     Write output to list
531
      if ( outmode.eq.-1 ) then
539
      if ( outmode.eq.-1 ) then
Line 614... Line 622...
614
      endif 
622
      endif 
615
 
623
 
616
c     Transform coordinates
624
c     Transform coordinates
617
      do i=1,npoints
625
      do i=1,npoints
618
 
626
 
-
 
627
          if ( ( abs(xx0(i)-mdv).gt.eps ).and.
-
 
628
     >         ( abs(yy0(i)-mdv).gt.eps )  )
-
 
629
     >    then
619
          rid=(xx0(i)-xmin)/dx+1.
630
              rid=(xx0(i)-xmin)/dx+1.
620
          rjd=(yy0(i)-ymin)/dy+1.
631
              rjd=(yy0(i)-ymin)/dy+1.
621
 
-
 
622
          xx0(i) = int_index3(lon,nx,ny,1,rid,rjd,1.,mdv)
632
              xx0(i) = int_index3(lon,nx,ny,1,rid,rjd,1.,mdv)
623
          yy0(i) = int_index3(lat,nx,ny,1,rid,rjd,1.,mdv)
633
              yy0(i) = int_index3(lat,nx,ny,1,rid,rjd,1.,mdv)
-
 
634
          else
-
 
635
              xx0(i) = mdv
-
 
636
              yy0(i) = mdv
-
 
637
          endif
624
 
638
 
625
      enddo
639
      enddo
626
 
640
 
627
c     Write output to list
641
c     Write output to list
628
      if ( outmode.eq.-1 ) then
642
      if ( outmode.eq.-1 ) then
Line 823... Line 837...
823
      endif
837
      endif
824
 
838
 
825
c     Transform coordinates
839
c     Transform coordinates
826
      do i=1,npoints
840
      do i=1,npoints
827
 
841
 
-
 
842
         if ( ( abs(xx0(i)-mdv).gt.eps ).and.
-
 
843
     >        ( abs(yy0(i)-mdv).gt.eps ).and.
-
 
844
     >        ( abs(zz0(i)-mdv).gt.eps ) )
-
 
845
     >    then
828
          call get_index3 (rid,rjd,rkd,xx0(i),yy0(i),-zz0(i),1,
846
             call get_index3 (rid,rjd,rkd,xx0(i),yy0(i),-zz0(i),1,
829
     >                     -p3,-ps,nx,ny,nz,xmin,ymin,dx,dy)
847
     >                        -p3,-ps,nx,ny,nz,xmin,ymin,dx,dy)
830
 
848
 
831
          zz0(i) = int_index3(z3,nx,ny,nz,rid,rjd,rkd,mdv)
849
             zz0(i) = int_index3(z3,nx,ny,nz,rid,rjd,rkd,mdv)
-
 
850
          else
-
 
851
             xx0(i) = mdv
-
 
852
             yy0(i) = mdv
-
 
853
             zz0(i) = mdv
-
 
854
          endif
832
 
855
 
833
      enddo
856
      enddo
834
 
857
 
835
c     Write output to list
858
c     Write output to list
836
      if ( outmode.eq.-1 ) then
859
      if ( outmode.eq.-1 ) then
Line 936... Line 959...
936
      endif
959
      endif
937
 
960
 
938
c     Transform coordinates
961
c     Transform coordinates
939
      do i=1,npoints
962
      do i=1,npoints
940
 
963
 
-
 
964
         if ( ( abs(xx0(i)-mdv).gt.eps ).and.
-
 
965
     >        ( abs(yy0(i)-mdv).gt.eps ).and.
-
 
966
     >        ( abs(zz0(i)-mdv).gt.eps ) )
-
 
967
     >    then
941
          call get_index3 (rid,rjd,rkd,xx0(i),yy0(i),zz0(i),1,
968
               call get_index3 (rid,rjd,rkd,xx0(i),yy0(i),zz0(i),1,
942
     >                     z3,zb,nx,ny,nz,xmin,ymin,dx,dy)
969
     >                          z3,zb,nx,ny,nz,xmin,ymin,dx,dy)
943
 
970
 
944
          zz0(i) = int_index3(p3,nx,ny,nz,rid,rjd,rkd,mdv)
971
               zz0(i) = int_index3(p3,nx,ny,nz,rid,rjd,rkd,mdv)
-
 
972
          else
-
 
973
               xx0(i) = mdv
-
 
974
               yy0(i) = mdv
-
 
975
               zz0(i) = mdv
-
 
976
          endif
945
 
977
 
946
      enddo
978
      enddo
947
 
979
 
948
c     Write output to list
980
c     Write output to list
949
      if ( outmode.eq.-1 ) then
981
      if ( outmode.eq.-1 ) then