module grid_module use input_module,only: | dlon_in,dlon_out, dlat_in,dlat_out, ! horizontal resolutions | vertical_in, vertical_out, ! vertical resolutions | interp_horizontal, interp_vertical implicit none type grid_type ! ! Geographic grid: integer :: nlat,nlon,nlev,nilev real :: glat1,glon1,zlev1,zilev1 real :: dlat,dlon,dlev real,pointer :: glat(:),glon(:), | zlev(:), ! midpoint level coord array (nlev) | zilev(:) ! interface level coord array (nilev) ! ! Magnetic grid: integer :: mlat,mlon,mlev,imlev real :: mlat1,mlon1,zmlev1,zimlev1,dmlev real,pointer :: gmlat(:),gmlon(:), | zmlev(:), ! midpoint level coord array (mlev) | zimlev(:) ! interface level coord array (imlev) end type grid_type type(grid_type) :: oldgrid,newgrid contains !----------------------------------------------------------------------- subroutine set_newgrid ! ! Set new grid (old grid is read from input history): ! integer :: ier,i,j,k ! ! Geographic grid: ! integer :: nlat,nlon,nlev,nilev ! real :: glat1,glon1,zlev1,zilev1 ! real :: dlat,dlon,dlev ! real,pointer :: glat(:),glon(:), ! | zlev(:), ! midpoint level coord array (nlev) ! | zilev(:) ! interface level coord array (nilev) ! ! New latitude grid: if (dlat_in /= dlat_out) then newgrid%dlat = dlat_out newgrid%glat1 = -90.+0.5*newgrid%dlat newgrid%nlat = int(abs(newgrid%glat1)*2./newgrid%dlat)+1 else newgrid%dlat = oldgrid%dlat newgrid%glat1 = oldgrid%glat1 newgrid%nlat = oldgrid%nlat endif ! if (associated(newgrid%glat)) deallocate(newgrid%glat) allocate(newgrid%glat(newgrid%nlat),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%glat:', | ' newgrid%nlat=',i3)") newgrid%nlat endif do j=1,newgrid%nlat newgrid%glat(j) = newgrid%glat1+(j-1)*newgrid%dlat enddo ! ! New longitude grid: if (dlon_in /= dlon_out) then newgrid%dlon = dlon_out newgrid%glon1 = -180. newgrid%nlon = 360./newgrid%dlon else newgrid%dlon = oldgrid%dlon newgrid%glon1 = oldgrid%glon1 newgrid%nlon = oldgrid%nlon endif if (associated(newgrid%glon)) deallocate(newgrid%glon) allocate(newgrid%glon(newgrid%nlon),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%glon:', | ' newgrid%nlon=',i3)") newgrid%nlon endif do i=1,newgrid%nlon newgrid%glon(i) = newgrid%glon1+(i-1)*newgrid%dlon enddo ! ! New levels grid: ! Midpoints levels dimensions, delta, and first: if (vertical_in /= vertical_out) then newgrid%dlev = vertical_out newgrid%zlev1 = oldgrid%zlev1-0.5*newgrid%dlev newgrid%nlev = int(((oldgrid%zlev(oldgrid%nlev)+.001)- | oldgrid%zlev1)/newgrid%dlev)+1 ! ! Interface levels dimensions, delta, and first: newgrid%zilev1 = oldgrid%zilev1 newgrid%nilev = int(((oldgrid%zilev(oldgrid%nilev)+.001)- | oldgrid%zilev1)/newgrid%dlev)+1 ! write(6,"('set_newgrid nilev: old nilev=',i3,' old zilev(nilev)=', ! | f8.2,' old zilev1=',f8.2,' new dlev=',f8.3,' new nilev=',i3)") ! | oldgrid%nilev,oldgrid%zilev(oldgrid%nilev),oldgrid%zilev1, ! | newgrid%dlev,newgrid%nilev else newgrid%dlev = oldgrid%dlev newgrid%zlev1 = oldgrid%zlev1 newgrid%zilev1 = oldgrid%zilev1 newgrid%nlev = oldgrid%nlev newgrid%nilev = oldgrid%nilev endif ! ! Midpoints levels coord arrays: if (associated(newgrid%zlev)) deallocate(newgrid%zlev) allocate(newgrid%zlev(newgrid%nlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%zlev:', | ' newgrid%nlev=',i3)") newgrid%nlev endif do k=1,newgrid%nlev newgrid%zlev(k) = newgrid%zlev1+(k-1)*newgrid%dlev enddo ! ! Interface levels coord arrays: if (associated(newgrid%zilev)) deallocate(newgrid%zilev) allocate(newgrid%zilev(newgrid%nilev),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%zilev:', | ' newgrid%nilev=',i3)") newgrid%nilev endif do k=1,newgrid%nilev newgrid%zilev(k) = newgrid%zilev1+(k-1)*newgrid%dlev enddo ! ! Magnetic grid: ! ! The vertical mag grid may be interpolated from 0.5 to 0.25, ! but the horizontal magnetic grid is unchanged from old history) ! ! integer :: mlat,mlon,mlev,imlev ! real :: mlat1,mlon1,zmlev1,zimlev1 ! real,pointer :: gmlat(:),gmlon(:), ! | zmlev(:), ! midpoint level coord array (mlev) ! | zimlev(:) ! interface level coord array (imlev) newgrid%mlat = oldgrid%mlat if (associated(newgrid%gmlat)) deallocate(newgrid%gmlat) allocate(newgrid%gmlat(newgrid%mlat),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%gmlat:', | ' newgrid%mlat=',i3)") newgrid%mlat endif newgrid%gmlat = oldgrid%gmlat newgrid%mlon = oldgrid%mlon if (associated(newgrid%gmlon)) deallocate(newgrid%gmlon) allocate(newgrid%gmlon(newgrid%mlon),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%gmlon:', | ' newgrid%mlon=',i3)") newgrid%mlon endif newgrid%gmlon = oldgrid%gmlon ! ! New mag levels grid: ! Midpoint mag levels dimensions, delta, and first: if (vertical_in /= vertical_out) then newgrid%dmlev = vertical_out newgrid%zmlev1 = oldgrid%zmlev1-0.5*newgrid%dmlev newgrid%mlev = int(((oldgrid%zmlev(oldgrid%mlev)+.001)- | oldgrid%zmlev1)/newgrid%dmlev)+1 ! ! Interface mag levels dimensions, delta, and first: newgrid%zimlev1 = oldgrid%zimlev1 newgrid%imlev = int(((oldgrid%zimlev(oldgrid%imlev)+.001)- | oldgrid%zimlev1)/newgrid%dmlev)+1 else ! ! 8/4/10 bugfix: Retain mag grid -- changed dlev to dmlev, etc. ! Thanks to Nicholas Pedatella newgrid%dmlev = oldgrid%dmlev newgrid%zmlev1 = oldgrid%zmlev1 newgrid%zimlev1 = oldgrid%zimlev1 newgrid%mlev = oldgrid%mlev newgrid%imlev = oldgrid%imlev endif ! ! Midpoints mag levels coord arrays: if (associated(newgrid%zmlev)) deallocate(newgrid%zmlev) allocate(newgrid%zmlev(newgrid%mlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%zmlev:', | ' newgrid%mlev=',i3)") newgrid%mlev endif do k=1,newgrid%mlev newgrid%zmlev(k) = newgrid%zmlev1+(k-1)*newgrid%dmlev enddo ! ! Interface mag levels coord arrays: if (associated(newgrid%zimlev)) deallocate(newgrid%zimlev) allocate(newgrid%zimlev(newgrid%imlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_newgrid: error allocating newgrid%zimlev:', | ' newgrid%imlev=',i3)") newgrid%imlev endif do k=1,newgrid%imlev newgrid%zimlev(k) = newgrid%zimlev1+(k-1)*newgrid%dmlev enddo ! ! Report to stdout: ! write(6,"(/,'set_newgrid: nlat=',i3,' dlat=',f6.2,' glat=',/, ! | (9f8.2))") newgrid%nlat,newgrid%dlat,newgrid%glat ! write(6,"(/,'set_newgrid: nlon=',i3,' dlon=',f6.2,' glon=',/, ! | (9f8.2))") newgrid%nlon,newgrid%dlon,newgrid%glon ! write(6,"(/,'set_newgrid: nlev=',i3,' dlev=',f6.2,' zlev=',/, ! | (9f8.3))") newgrid%nlev,newgrid%dlev,newgrid%zlev ! write(6,"(/,'set_newgrid: nilev=',i3,' dlev=',f6.2,' zilev=',/, ! | (9f8.3))") newgrid%nilev,newgrid%dlev,newgrid%zilev ! write(6,"(/,'set_newgrid: mlat=',i3,' gmlat=',/,(9f8.2))") ! | newgrid%mlat,newgrid%gmlat ! write(6,"(/,'set_newgrid: mlon=',i3,' gmlon=',/,(9f8.2))") ! | newgrid%mlon,newgrid%gmlon ! write(6,"(/,'set_newgrid: mlev=',i3,' gmlev=',/,(9f8.3))") ! | newgrid%mlev,newgrid%zmlev ! write(6,"(/,'set_newgrid: imlev=',i3,' gimlev=',/,(9f8.3))") ! | newgrid%imlev,newgrid%zimlev ! write(6,"(' ')") end subroutine set_newgrid !----------------------------------------------------------------------- subroutine set_oldgrid(h) use readhist_module,only: history type(history),intent(in) :: h ! integer :: ier ! ! Latitude: oldgrid%nlat = h%nlat allocate(oldgrid%glat(oldgrid%nlat),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%glat:', | ' oldgrid%nlat=',i3)") oldgrid%nlat endif oldgrid%glat(:) = h%glat(:) oldgrid%glat1 = h%glat(1) oldgrid%dlat = oldgrid%glat(2)-oldgrid%glat(1) ! assume regular spacing ! ! Longitude: oldgrid%nlon = h%nlon allocate(oldgrid%glon(oldgrid%nlon),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%glon:', | ' oldgrid%nlon=',i3)") oldgrid%nlon endif oldgrid%glon(:) = h%glon(:) oldgrid%glon1 = h%glon(1) oldgrid%dlon = oldgrid%glon(2)-oldgrid%glon(1) ! assume regular spacing ! ! Midpoint Levels: oldgrid%nlev = h%nlev allocate(oldgrid%zlev(oldgrid%nlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%zlev:', | ' oldgrid%nlev=',i3)") oldgrid%nlev endif oldgrid%zlev(:) = h%lev(:) oldgrid%zlev1 = h%lev(1) oldgrid%dlev = oldgrid%zlev(2)-oldgrid%zlev(1) ! assume regular spacing ! ! Interface Levels: oldgrid%nilev = h%nilev allocate(oldgrid%zilev(oldgrid%nilev),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%zilev:', | ' oldgrid%nilev=',i3)") oldgrid%nilev endif oldgrid%zilev(:) = h%ilev(:) oldgrid%zilev1 = h%ilev(1) oldgrid%dlev = oldgrid%zilev(2)-oldgrid%zilev(1) ! assume regular spacing ! ! Magnetic grid: ! ! integer :: mlat,mlon,mlev,imlev ! real :: mlat1,mlon1,zmlev1,zimlev1 ! real,pointer :: gmlat(:),gmlon(:), ! | zmlev(:), ! midpoint level coord array (mlev) ! | zimlev(:) ! interface level coord array (imlev) ! ! ! Mag Latitude (same as input history): oldgrid%mlat = h%mlat allocate(oldgrid%gmlat(oldgrid%mlat),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%gmlat:', | ' oldgrid%mlat=',i3)") oldgrid%mlat endif oldgrid%gmlat(:) = h%gmlat(:) ! ! Mag Longitude (same as old grid): oldgrid%mlon = h%mlon allocate(oldgrid%gmlon(oldgrid%mlon),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%gmlon:', | ' oldgrid%mlon=',i3)") oldgrid%mlon endif oldgrid%gmlon(:) = h%gmlon(:) ! ! Mag Midpoint Levels: oldgrid%mlev = h%mlev allocate(oldgrid%zmlev(oldgrid%mlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%zmlev:', | ' oldgrid%mlev=',i3)") oldgrid%mlev endif oldgrid%zmlev(:) = h%gmlev(:) oldgrid%zmlev1 = h%gmlev(1) ! ! Mag Interface Levels: oldgrid%imlev = h%imlev allocate(oldgrid%zimlev(oldgrid%imlev),stat=ier) if (ier /= 0) then write(6,"('>>> set_oldgrid: error allocating oldgrid%zimlev:', | ' oldgrid%imlev=',i3)") oldgrid%imlev endif oldgrid%zimlev(:) = h%gimlev(:) oldgrid%zimlev1 = h%gimlev(1) oldgrid%dmlev = oldgrid%zmlev(2)-oldgrid%zmlev(1) ! assume regular spacing ! ! Report to stdout: ! write(6,"(/,'set_oldgrid: nlat=',i3,' dlat=',f6.2,' glat=',/, ! | (9f8.2))") oldgrid%nlat,oldgrid%dlat,oldgrid%glat ! write(6,"(/,'set_oldgrid: nlon=',i3,' dlon=',f6.2,' glon=',/, ! | (9f8.2))") oldgrid%nlon,oldgrid%dlon,oldgrid%glon ! write(6,"(/,'set_oldgrid: nlev=',i3,' dlev=',f6.2,' zlev=',/, ! | (9f8.2))") oldgrid%nlev,oldgrid%dlev,oldgrid%zlev ! write(6,"(/,'set_oldgrid: nilev=',i3,' dlev=',f6.2,' zilev=',/, ! | (9f8.2))") oldgrid%nilev,oldgrid%dlev,oldgrid%zilev ! write(6,"(/,'set_oldgrid: mlat=',i3,' gmlat=',/,(9f8.2))") ! | oldgrid%mlat,oldgrid%gmlat ! write(6,"(/,'set_oldgrid: mlon=',i3,' gmlon=',/,(9f8.2))") ! | oldgrid%mlon,oldgrid%gmlon ! write(6,"(/,'set_oldgrid: mlev=',i3,' zmlev=',/,(9f8.3))") ! | oldgrid%mlev,oldgrid%zmlev ! write(6,"(/,'set_oldgrid: imlev=',i3,' zimlev=',/,(9f8.3))") ! | oldgrid%imlev,oldgrid%zimlev ! write(6,"(' ')") end subroutine set_oldgrid end module grid_module