C********************************************************************* C C PROGRAM xzcont_uwtq_cl C C x-z平面等値線図作成プログラム(時刻, 条件指定, カラー) C C 履歴 1997/07/24 中島健介 C 1997/07/25 小高正嗣 C 1997/07/26 小高正嗣 C 1997/07/26 小高正嗣 C 1998/09/15 小高正嗣 C 1998/10/02 小高正嗣 C 1999/01/17 小高正嗣 C 1999/01/18 小高正嗣 C 1999/01/25 小高正嗣 C 1999/02/04 小高正嗣 C C********************************************************************* program xzcont_uwt_cl *----------------------------------------------------- implicit real*8 ( a-h, o-z ) *----------------------------------------------------- #include "grid_size_M.f" parameter ( nbr3 = 3*((nx+6)*(nz+6)/16+1) ) *----------------------------------------------------- integer ibr(nbr3) *----------------------------------------------------- parameter ( day = 86400.0 ) parameter ( hour = 3600.0 ) parameter ( PI = 3.141592D0 ) real cdp_uw,cdp_pt,cdp_qv,cdp_ptd ! カラートーン刻み real tlev1,tlev2 ! カラーレベル integer ipat ! トーン番号 integer cmax parameter ( cmax = 60 ) ! カラー数 real wmax, wmin, umax, umin, ptmax, ptmin, qvmax, qvmin real ptdmax, ptdmin parameter ( wmax = 25., wmin = -25., umax = 25., umin = -25. ) parameter ( ptmax = 280., ptmin = 180. ) parameter ( ptdmax = 4.0, ptdmin = -4.0 ) parameter ( qvmax = 10.e-6, qvmin = 1.e-6 ) real*4 pti_u(2,cmax+1), pti_w(2,cmax+1) real*4 pti_pt(2,cmax+1), pti_ptd(2,cmax+1) real*4 pti_qv(2,cmax+1) *----------------------------------------------------- real*8 \ zgrid ( -2:nzb ), \ temp0 ( -2:nzb ), \ tpot0 ( -2:nzb ), \ dens0 ( -2:nzb ), \ ppai0 ( -2:nzb ), \ pres0 ( -2:nzb ), \ rvap0 ( -2:nzb ), \ qvap0 ( -2:nzb ) real*8 \ u ( -2:nxb, -2:nzb ), \ w ( -2:nxb, -2:nzb ), \ tpot ( -2:nxb, -2:nzb ), \ qvap ( -2:nxb, -2:nzb ) real*4 \ z0p ( 0:nz ), \ p0p ( 0:nz ), \ d0p ( 0:nz ), \ temp0p ( 0:nz ), \ tpot0p ( 0:nz ), \ pai0p ( 0:nz ), \ rvap0p ( 0:nz ), \ qvap0p ( 0:nz ) real*4 \ up ( nx+1,nz+1 ), \ wp ( nx+1,nz+1 ), \ tpotp ( nx+1, nz+1 ), \ tpall ( nx+1, nz+1 ), \ tpave ( nz+1 ), \ tpdev ( nx+1, nz+1 ), \ qvapp ( nx+1, nz+1 ) character*2 mint character*2 lt character*2 date character*17 header character*100 dumy character CSGI*1 *----------------------------------------------------- open (11,file='exparam') read (11,*) irun read (11,*) time0, nloop1, nloop2, dtime, ntmoni, ntkubu, times read (11,*) ylat, als call FOPEN ( irun+1 ) do iz = 1,9 read(16,'(a100)') dumy end do do iz = -2, nzb read(16,*) ii, zgrid(iz), pres0(iz), dens0(iz), temp0(iz), \ tpot0(iz), ppai0(iz), qvap0(iz), rvap0(iz) end do do iz = 0, nz z0p(iz) = zgrid(iz) p0p(iz) = pres0(iz) d0p(iz) = dens0(iz) temp0p(iz)= temp0(iz) tpot0p(iz)= tpot0(iz) pai0p(iz) = ppai0(iz) rvap0p(iz)= rvap0(iz) qvap0p(iz)= qvap0(iz) end do write(6,*) 'data number = ', nloop1 write(6,*) 'INPUT DATA BEGINING, END, & INTERVAL' write(6,*) 'itb ite interval ?' read (*,*) itb, ite, interv write(6,*) 'SELECT TEMP. MODE (I)' write(6,*) '1; RAW PTEMP 2; DEVIATION ?' read (5,*) IDN call sgpwsn read(*,*) iws c CALL SWISET('IWIDTH', 900) c CALL SWISET('IHEIGHT', 400) c call swlset( 'LDUMP', .TRUE. ) call gropn(iws) call sldiv( 'T', 2, 2 ) c call gropn(-iws) c call sldiv( 'T', 1, 3 ) call swlset( 'LALT', .TRUE. ) call sglset( 'LFULL', .TRUE. ) call sglset( 'LCORNER', .FALSE. ) call uzfact(0.6) call udlset( 'LMSG', .FALSE. ) call udiset( 'INDXMJ', 2 ) cdp_uw = (wmax - wmin)/cmax cdp_pt = (ptmax - ptmin)/cmax cdp_qv = (qvmax - qvmin)/cmax cdp_ptd = (ptdmax - ptdmin)/cmax do k=1,cmax+1 pti_w(1,k) = wmin + (k-1)*cdp_uw pti_w(2,k) = wmin + (k-1)*cdp_uw pti_pt(1,k) = ptmin + (k-1)*cdp_pt pti_pt(2,k) = ptmin + (k-1)*cdp_pt pti_qv(1,k) = qvmin + (k-1)*cdp_qv pti_qv(2,k) = qvmin + (k-1)*cdp_qv pti_ptd(1,k) = ptdmin + (k-1)*cdp_ptd pti_ptd(2,k) = ptdmin + (k-1)*cdp_ptd end do do it = 1, ite read(31) time, ((u(ix,iz),ix=-2,nxb), iz=-2,nzb) read(33) time, ((w(ix,iz),ix=-2,nxb), iz=-2,nzb) read(34) time, ((tpot(ix,iz),ix=-2,nxb), iz=-2,nzb) read(35) time, ((qvap(ix,iz),ix=-2,nxb), iz=-2,nzb) write(0,*) 'reading time = ', time if ( it .ge. itb .and. it .le. ite ) then if ( mod( (it-itb), interv ) .eq. 0 ) then write(0,*) 'writing time = ', time c*** 時刻変換: 秒から LT へ************************************ time = time + times write(date,'(I2.2)') INT( time / day ) + 1 time = ( time - INT( time / day ) * day ) / hour if ( INT(time) .ge. 24 ) then time = time - 24 end if write(lt, '(I2.2)') INT(time) write(mint, '(I2.2)') INT( MOD(time,1)*60.0 ) header = 'day=' // date // ' ' // 'LT=' // lt // ':' // mint c************************************************************** c************************************************************** do iz = 1,nz+1 do ix = 1, nx+1 up(ix,iz) = u(ix-1,iz-1) c wp(ix,iz) = w(ix-1,iz-1) wp(ix,iz) = ( w(ix-2,iz-1) + \ 2.0*w(ix-1,iz-1) + w(ix,iz-1))/4.0 qvapp(ix,iz) = qvap(ix-1,iz-1) c tpotp(ix,iz) = tpot(ix-1,iz-1) tpotp(ix,iz) = ( tpot(ix-2,iz-1) + \ 2.0*tpot(ix-1,iz-1) + tpot(ix,iz-1))/4.0 c tpall(ix,iz) = tpot(ix-1,iz-1) + tpot0(iz-1) tpall(ix,iz) = tpotp(ix,iz) + tpot0(iz-1) end do end do do iz = 1, nz+1 tpave(iz) = 0. do ix = 1,nx+1 tpave(iz) = tpave(iz) + tpall(ix,iz) /(nx+1) end do end do do iz = 1,nz+1 do ix = 1, nx+1 tpdev(ix,iz) = tpall(ix,iz) - tpave(iz) end do end do if ( IDN .eq. 1 ) then C********************************************************************* C pt C********************************************************************* do i = 1,cmax TLEV2 = ptmin + i*cdp_pt TLEV1 = TLEV2 - cdp_pt ipat = (30+(i-1))*1000 + 999 call uestlv(TLEV1, TLEV2, IPAT) end do call udiset( 'ICYCLE', 5 ) call grfrm call grfig call grswnd ( 0.0, 51.2, 0.0, 20.0 ) call grsvpt ( 0.1, 0.9, 0.1, 0.42 ) c call grsvpt ( 0.1, 0.9, 0.15, 0.306 ) call grstrn (1) call uspfit call grstrf call udpset ( 'rsizet', 0.01 ) call ussttl ( CSGI(135)//'[K]', 'X[km]', '', 'Z[km]' ) call uetonf ( tpall, nx+1, nx+1, nz+1 ) call usdaxs call uxsttl ( 'T', header, 1. ) C********************************************************************* C トーンバー C********************************************************************* call grfig call grswnd( 0.0, 1.0, ptmin, ptmax ) call grsvpt ( 0.95, 0.98, 0.1, 0.42) c call grsvpt( 0.95, 0.98, 0.15, 0.306 ) call grstrn( 1 ) call grstrf call uetone( pti_pt, 2, 2, cmax+1 ) call slpvpr( 3 ) call uyaxdv( 'L', 2., 10. ) call uzinit call ueitlv else C********************************************************************* C ptemp_dev C********************************************************************* call udiset( 'ICYCLE', 10 ) TLEV2 = ptdmin + cdp_ptd TLEV1 = -10.0 call uestlv(TLEV1, TLEV2, 30999) c do i = 1,cmax do i = 2,cmax-1 TLEV2 = ptdmin + i*cdp_ptd TLEV1 = TLEV2 - cdp_ptd ipat = (30+(i-1))*1000 + 999 call uestlv(TLEV1, TLEV2, IPAT) end do TLEV2 = 10.0 TLEV1 = ptdmax - cdp_ptd call uestlv(TLEV1, TLEV2, 89999) call grfrm call grfig call grswnd ( 0.0, 51.2, 0.0, 20.0 ) call grsvpt ( 0.1, 0.9, 0.1, 0.42 ) c call grsvpt ( 0.1, 0.9, 0.15, 0.306 ) call grstrn (1) call uspfit call grstrf call udpset ( 'rsizet', 0.01 ) call ussttl ( CSGI(131) \ //CSGI(135)//'[K]', 'X[km]', '', 'Z[km]' ) call uetonf ( tpdev, nx+1, nx+1, nz+1 ) call usdaxs call uxsttl ( 'T', header, 1. ) c call udgcla ( -10.0, 10.0, 2.0 ) c call uddclv ( 2.0 ) c call uddclv ( -2.0 ) c call udsclv ( 0.0, 1, 1, '', '0.0' ) c call udcntz ( tpdev, nx+1, nx+1, nz+1, ibr, nbr3) C********************************************************************* C トーンバー C********************************************************************* call grfig call grswnd( 0.0, 1.0, ptdmin, ptdmax ) call grsvpt ( 0.95, 0.98, 0.1, 0.42) c call grsvpt( 0.95, 0.98, 0.15, 0.306 ) call grstrn( 1 ) call grstrf call uetone( pti_ptd, 2, 2, cmax+1 ) call slpvpr( 3 ) call uyaxdv( 'L', 0.2, 1.0 ) call uzinit call ueitlv end if * * C********************************************************************* C w C********************************************************************* call udiset( 'ICYCLE', 4 ) TLEV2 = wmin TLEV1 = -50.0 call uestlv(TLEV1, TLEV2, 30999) cdp_uw = (wmax - wmin)/(cmax-2) do i = 1,cmax-2 TLEV2 = wmin + i*cdp_uw TLEV1 = TLEV2 - cdp_uw ipat = (30+i)*1000 + 999 call uestlv(TLEV1, TLEV2, IPAT) end do TLEV2 = 50.0 TLEV1 = wmax call uestlv(TLEV1, TLEV2, 89999) call grfrm call grfig call grswnd ( 0.0, 51.2, 0.0, 20.0 ) call grsvpt ( 0.1, 0.9, 0.1, 0.42 ) c call grsvpt ( 0.1, 0.9, 0.15, 0.306 ) call grstrn (1) call uspfit call grstrf call udpset ( 'rsizet', 0.01 ) call ussttl ( 'W[m/sec]', 'X[km]', '', 'Z[km]' ) call uetonf ( wp, nx+1, nx+1, nz+1 ) call usdaxs call uxsttl ( 'T', header, 1. ) call udgcla ( -50.0, 50.0, 5.0 ) call udsclv ( 0.0, 1, 1, '', '0.0' ) call udcntz ( wp, nx+1, nx+1, nz+1, ibr, nbr3) C********************************************************************* C トーンバー C********************************************************************* call grfig call grswnd( 0.0, 1.0, wmin, wmax ) call grsvpt ( 0.95, 0.98, 0.1, 0.42 ) c call grsvpt( 0.95, 0.98, 0.15, 0.306 ) call grstrn( 1 ) call grstrf call uetone( pti_w, 2, 2, cmax+1 ) call slpvpr( 3 ) call uyaxdv( 'L', 5., 10. ) call uzinit C********************************************************************* C u C********************************************************************* call grfrm call grfig call grswnd ( 0.0, 51.2, 0.0, 20.0 ) call grsvpt ( 0.1, 0.9, 0.1, 0.42 ) c call grsvpt ( 0.1, 0.9, 0.15, 0.306 ) call grstrn (1) call uspfit call grstrf call udpset ( 'rsizet', 0.01 ) call ussttl ( 'U[m/sec]', 'X[km]', '', 'Z[km]' ) call uetonf ( up, nx+1, nx+1, nz+1 ) call usdaxs call uxsttl ( 'T', header, 1. ) call udgcla ( -50.0, 50.0, 5.0 ) call udsclv ( 0.0, 1, 1, '', '0.0' ) call udcntz ( up, nx+1, nx+1, nz+1, ibr, nbr3) C********************************************************************* C トーンバー C********************************************************************* call grfig call grswnd( 0.0, 1.0, umin, umax ) call grsvpt ( 0.95, 0.98, 0.1, 0.42) c call grsvpt( 0.95, 0.98, 0.15, 0.306 ) call grstrn( 1 ) call grstrf call uetone( pti_w, 2, 2, cmax+1 ) call slpvpr( 3 ) call uyaxdv( 'L', 5., 10. ) call uzinit call ueitlv * * end if end if end do call grcls call FCLOSE end