program permutation implicit none integer::s,i,j,h,imn,imx,jmn,jmx,syear,eyear,iyear,smonth,emonth,imonth,ice1,jce1,xst,xen,yst,yen,month1,n,nto,nto2,nmx,nt,nnt,out,nout,nto3,it,itmax,tt,kk,m,nn,num integer::xst2,xen2,yst2,yen2,xst3,xen3,yst3,yen3,time,time2,time4,itime,t,tmx,tmn,ntmx,nntmx,trid,trid2,xst1,xen1,yst1,yen1,iit,nnew integer::kst, ken integer::hmx,hmn,hhmn,hhmx,hmx2 integer :: hlev integer::k,ememax,emem,emem1 parameter(imx=51,jmx=51,hmx=37,hmn=1,jmn=1,imn=1,tmn=1,hhmn=1,hhmx=40,hmx2=27,nntmx=10000,itmax=5000,nmx=100000,ememax=3000) integer,dimension(itmax)::nn2 integer, dimension(imx, jmx) :: nu, nv, nvor, ntmp, nspfh, nhgt, naetoke !integer,dimension(imx,jmx,hmx)::nv,nvvel,ntmp,npv,nhgt,nvor !integer,dimension(imx,jmx,hmx2)::nsp real::lat1,lon1,jlat,jlon,pi,es,ssp,slpup,slplow,lon2,lat2,minslp,hoge,varvn,varsrehn,s12,varcapen,varpvn,varspn,varhgtn,varvveln,vartmpn,rmin,rmin2 real,dimension(itmax) :: r,rr real,dimension(nntmx) :: trid3,time3 real,dimension(jmx) :: co,lat real,dimension(imx,jmx,hmx) :: spfh, u, v, vor, tmp, hgt real,dimension(imx,jmx) :: epe, depe, epeflux, baroc, diab, corre, aetoke real, dimension(imx, jmx, itmax) :: coepe, coepe2, coepe3 real, dimension(imx, jmx, itmax) :: codepe, codepe2, codepe3 real, dimension(imx, jmx, itmax) :: coepeflux, coepeflux2, coepeflux3 real, dimension(imx, jmx, itmax) :: cobaroc, cobaroc2, cobaroc3 real, dimension(imx, jmx, itmax) :: codiab, codiab2, codiab3 real, dimension(imx, jmx, itmax) :: cocorre, cocorre2, cocorre3 real, dimension(imx, jmx, itmax) :: coaetoke, coaetoke2, coaetoke3 real, dimension(imx, jmx) :: coepe4, codepe4, coepeflux4, cobaroc4 real, dimension(imx, jmx) :: codiab4, cocorre4, coaetoke4 real,dimension(jmx,itmax)::sco3,sco4,sco5 real, dimension(imx, jmx) :: aveu1, aveu2, du, dcou, pu real, dimension(imx, jmx) :: avev1, avev2, dv, dcov, pv real, dimension(imx, jmx) :: avevor1, avevor2, dvor, dcovor, pvor real, dimension(imx, jmx) :: avetmp1, avetmp2, dtmp, dcotmp, ptmp real, dimension(imx, jmx) :: avespfh1, avespfh2, dspfh, dcospfh, pspfh real, dimension(imx, jmx) :: avehgt1, avehgt2, dhgt, dcohgt, phgt real, dimension(imx, jmx) :: aveaetoke1, aveaetoke2, daetoke, dcoaetoke, paetoke real, dimension(imx, jmx, ememax) :: epe2d, depe2d, epeflux2d, baroc2d, diab2d, corre2d, aetoke2d real, dimension(imx, jmx, ememax) :: su1, su2, sv1, sv2, svor1, svor2, stmp1, stmp2 real, dimension(imx, jmx, ememax) :: svvel1, svvel2, sspfh1, sspfh2, spv1, spv2 ,shgt1,shgt2 real, dimension(imx, jmx, ememax) :: saetoke1, saetoke2 real::snum,snum1,snum2,snum3,tnum,tnum2,ss,snum11,snum12,snum13,snum21,snum22,snum23 real::p(hmx)=(/1000., 975., 950., 925., 900., 875., 850., 825., 800., 775., 750., 700., & 650., 600., 550., 500., 450., 400., 350., 300., 250., 225., 200., 175., 150., 125., 100., & 70., 50., 30., 20., 10.,7., 5., 3., 2., 1./) character::kt1*1,mmonth*2,year*4,SEA*3,region*3,kt*2,CAT*2,CAT2*2 character :: jetpos1*8, jetpos2*8 logical::fast !itime = 13 !CAT='PO' !CAT2='OJ' SEA='SON' region='JPN' su1 = 0 sv1 = 0 svor1 = 0 stmp1 = 0 sspfh1 = 0 shgt1 = 0 saetoke1 = 0 su2 = 0 sv2 = 0 svor2 = 0 stmp2 = 0 sspfh2 = 0 shgt2 = 0 saetoke2 = 0 snum11 = 0 snum12 = 0 snum13 = 0 snum21 = 0 snum22 = 0 snum23 = 0 coepe = 0 coepe2 = 0 codepe = 0 codepe2 = 0 coepeflux = 0 coepeflux2 = 0 cobaroc = 0 cobaroc2 = 0 codiab = 0 codiab2 = 0 cocorre = 0 cocorre2 = 0 coaetoke = 0 coaetoke2 = 0 jetpos1 = 'Njetconf' jetpos2 = 'Njetdiff' !write(6,*) 'read namelist' !namelist/namdim/CAT,CAT2,itime !write(6,*) 'category' !write(6,*) CAT,' ',CAT2,' ',itime !read(5,namdim) CAT='PO' itime = 9 kst = 3 ken = 3 hlev = 14 it = itime if(it.le.9) then write(kt1,'(i1)') it kt='0'//kt1 else write(kt,'(i2)') it endif open(69,file='pvalue_for_permutation_ecjet_epeterms_rotate_'//CAT//'_'//kt//'_'//jetpos1//'_'//jetpos2//'.data',status='unknown',form='unformatted') emem=1 do k=kst,ken ss = 1 if(k.eq.1) then SEA = 'MAM' elseif(k.eq.2) then SEA = 'SON' elseif(k.eq.3) then SEA = 'DJF' endif write(6, *) SEA open(68,file='./'//CAT//'/emem_'//kt//'_energy_notmean_1000_700_DRATE05_100_'//CAT//'_'//SEA//'_'//jetpos1//'_160E_190E.data',status='unknown',form='unformatted') !open(68,file=''//CAT//'/emem_exspfront600_ktvar_DRATE05_100_'//CAT//'_'//SEA//'_'//kt//'.data',status='unknown',form='unformatted') kk=68 do it=1,ememax read(kk,end=150) !((epe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((epe(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((epeflux(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((baroc(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((diab(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((corre(i,j),i=imn,imx),j=jmn,jmx) h = hlev do j = jmn, jmx do i = imn, imx epe2d(i,j,emem) = epe(i, j) depe2d(i,j,emem) = depe(i, j) epeflux2d(i,j,emem) = epeflux(i,j) baroc2d(i,j,emem) = baroc(i,j) diab2d(i,j,emem) = diab(i,j) corre2d(i,j,emem) = corre(i,j) aetoke2d(i,j,emem) = aetoke(i,j) enddo enddo !read(kk,end=150) ((sreh(i,j,emem),i=imn,imx),j=jmn,jmx) !read(kk) ((cold(i,j,emem),i=imn,imx),j=jmn,jmx) !read(kk) ((warm(i,j,emem),i=imn,imx),j=jmn,jmx) emem = emem + 1 ss = ss + 1 enddo 150 continue if(k.eq.1) then snum11 = ss endif if(k.eq.2) then snum12 = ss endif if(k.eq.3) then snum13 = ss endif enddo !kloop emem1 = emem write(6,*) 'emem1', emem1 do k=kst,ken ss = 1 if(k.eq.1) then SEA = 'MAM' elseif(k.eq.2) then SEA = 'SON' elseif(k.eq.3) then SEA = 'DJF' endif !open(68,file=''//CAT//'/emem_exnspfront600_ktvar_DRATE05_100_'//CAT//'_'//SEA//'_'//kt//'.data',status='unknown',form='unformatted') open(68,file='./'//CAT//'/emem_'//kt//'_energy_notmean_1000_700_DRATE05_100_'//CAT//'_'//SEA//'_'//jetpos2//'_160E_190E.data',status='unknown',form='unformatted') kk=68 do it=1,ememax read(kk,end=151) !((epe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) !((depe(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((epe(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((epeflux(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((baroc(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((diab(i,j),i=imn,imx),j=jmn,jmx) read(kk) ((corre(i,j),i=imn,imx),j=jmn,jmx) do j = jmn, jmx do i = imn, imx epe2d(i,j,emem) = epe(i, j) depe2d(i,j,emem) = depe(i, j) epeflux2d(i,j,emem) = epeflux(i,j) baroc2d(i,j,emem) = baroc(i,j) diab2d(i,j,emem) = diab(i,j) corre2d(i,j,emem) = corre(i,j) aetoke2d(i,j,emem) = aetoke(i,j) enddo enddo !read(kk,end=151) ((sreh(i,j,emem),i=imn,imx),j=jmn,jmx) !read(kk) ((cold(i,j,emem),i=imn,imx),j=jmn,jmx) !read(kk) ((warm(i,j,emem),i=imn,imx),j=jmn,jmx) emem = emem + 1 ss = ss + 1 enddo 151 continue if(k.eq.1) then snum21 = ss endif if(k.eq.2) then snum22 = ss endif if(k.eq.3) then snum23 = ss endif enddo !kloop write(6,*) 'emem', emem do s=kst,ken if(s.eq.1) then SEA='MAM' snum=snum11 elseif(s.eq.2) then SEA='SON' snum=snum12 else SEA='DJF' snum=snum13 endif open(57,file='compo_energy_notmean_1000_700_DRATE05_100_'//CAT//'_'//SEA//'_'//jetpos1//'_160E_190E.data',status='unknown',form='unformatted') tnum=snum11+snum12+snum13 write(6,*) snum11,snum12,snum13 do it=1,17 !AO read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) ((coepe4(i, j),i=imn,imx),j=imn,imx) read(57) ((coepeflux4(i, j),i=imn,imx),j=imn,imx) read(57) ((cobaroc4(i, j),i=imn,imx),j=imn,imx) read(57) ((codiab4(i, j),i=imn,imx),j=imn,imx) read(57) ((cocorre4(i, j),i=imn,imx),j=imn,imx) do j = jmn, jmx do i = imn, imx coepe3(i,j,it) = coepe4(i, j) codepe3(i,j,it) = codepe4(i, j) coepeflux3(i,j,it) = coepeflux4(i, j) cobaroc3(i,j,it) = cobaroc4(i, j) codiab3(i,j,it) = codiab4(i, j) cocorre3(i,j,it) = cocorre4(i, j) coaetoke3(i,j,it) = coaetoke4(i, j) enddo enddo !read(57) ((cosreh3(i,j,it),i=imn,imx),j=jmn,jmx) !read(57) ((cocold3(i,j,it),i=imn,imx),j=jmn,jmx) !read(57) ((cowarm3(i,j,it),i=imn,imx),j=jmn,jmx) do j=jmn,jmx do i=imn,imx coepe(i,j,it) = coepe(i,j,it) + coepe3(i,j,it)*snum/tnum codepe(i,j,it) = codepe(i,j,it) + codepe3(i,j,it)*snum/tnum coepeflux(i,j,it) = coepeflux(i,j,it) + coepeflux3(i,j,it)*snum/tnum cobaroc(i,j,it) = cobaroc(i,j,it) + cobaroc3(i,j,it)*snum/tnum codiab(i,j,it) = codiab(i,j,it) + codiab3(i,j,it)*snum/tnum cocorre(i,j,it) = cocorre(i,j,it) + cocorre3(i,j,it)*snum/tnum coaetoke(i,j,it) = coaetoke(i,j,it) + coaetoke3(i,j,it)*snum/tnum enddo enddo enddo if(s.eq.1) then SEA='MAM' snum=snum21 elseif(s.eq.2) then SEA='SON' snum=snum22 else SEA='DJF' snum=snum23 endif open(58,file='compo_energy_notmean_1000_700_DRATE05_100_'//CAT//'_'//SEA//'_'//jetpos2//'_160E_190E.data',status='unknown',form='unformatted') !write(6, *) CAT, SEA tnum2=snum21+snum22+snum23 write(6,*) snum21,snum22,snum23 do it=1, 17 !PO read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) !((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) ((coepe4(i, j),i=imn,imx),j=imn,imx) read(58) ((coepeflux4(i, j),i=imn,imx),j=imn,imx) read(58) ((cobaroc4(i, j),i=imn,imx),j=imn,imx) read(58) ((codiab4(i, j),i=imn,imx),j=imn,imx) read(58) ((cocorre4(i, j),i=imn,imx),j=imn,imx) do j = jmn, jmx do i = imn, imx coepe3(i,j,it) = coepe4(i, j) codepe3(i,j,it) = codepe4(i, j) coepeflux3(i,j,it) = coepeflux4(i, j) cobaroc3(i,j,it) = cobaroc4(i, j) codiab3(i,j,it) = codiab4(i, j) cocorre3(i,j,it) = cocorre4(i, j) coaetoke3(i,j,it) = coaetoke4(i, j) enddo enddo do j=jmn,jmx do i=imn,imx coepe2(i,j,it) = coepe2(i,j,it) + coepe3(i,j,it)*snum/tnum2 codepe2(i,j,it) = codepe2(i,j,it) + codepe3(i,j,it)*snum/tnum2 coepeflux2(i,j,it) = coepeflux2(i,j,it) + coepeflux3(i,j,it)*snum/tnum2 cobaroc2(i,j,it) = cobaroc2(i,j,it) + cobaroc3(i,j,it)*snum/tnum2 codiab2(i,j,it) = codiab2(i,j,it) + codiab3(i,j,it)*snum/tnum2 cocorre2(i,j,it) = cocorre2(i,j,it) + cocorre3(i,j,it)*snum/tnum2 coaetoke2(i,j,it) = coaetoke2(i,j,it) + coaetoke3(i,j,it)*snum/tnum2 enddo enddo enddo enddo !season loop it= itime do j=jmn,jmx do i=imn,imx dcou(i,j) = coepe(i,j,it) - coepe2(i,j,it) dcov(i,j) = codepe(i,j,it) - codepe2(i,j,it) dcovor(i,j) = coepeflux(i,j,it) - coepeflux2(i,j,it) dcotmp(i,j) = cobaroc(i,j,it) - cobaroc2(i,j,it) dcospfh(i,j) = codiab(i,j,it) - codiab2(i,j,it) dcohgt(i,j) = cocorre(i,j,it) - cocorre2(i,j,it) dcoaetoke(i,j) = coaetoke(i,j,it) - coaetoke2(i,j,it) enddo enddo nu = 0 nv = 0 nvor = 0 ntmp = 0 nspfh = 0 nhgt = 0 naetoke = 0 do nnt=1,nmx aveu1 = 0. avev1 = 0. avevor1 = 0. avetmp1 = 0. avespfh1 = 0. avehgt1 = 0. aveaetoke1 = 0. aveu2 = 0. avev2 = 0. avevor2 = 0. avetmp2 = 0. avespfh2 = 0. avehgt2 = 0. aveaetoke2 = 0. rmin2=0. do n=1,emem call random_seed call random_number(r(1:emem)) rmin=1. do nn=1,emem if(r(nn).lt.rmin.and.r(nn).gt.rmin2) then rmin=r(nn) nnew=nn endif enddo rr(n)=rmin nn2(n)=nnew rmin2=rmin enddo !do m=1,itmax !write(6,*) rr(m),nn2(m) !enddo do it=1,emem do j=jmn,jmx do i=imn,imx if(it.le.emem1) then su1(i,j,it) = epe2d(i, j, nn2(it)) sv1(i,j,it) = depe2d(i, j, nn2(it)) svor1(i,j,it) = epeflux2d(i, j, nn2(it)) stmp1(i,j,it) = baroc2d(i, j, nn2(it)) sspfh1(i,j,it) = diab2d(i, j, nn2(it)) shgt1(i,j,it) = corre2d(i, j, nn2(it)) saetoke1(i,j,it) = aetoke2d(i, j, nn2(it)) else su2(i,j,it - emem1) = epe2d(i, j, nn2(it)) sv2(i,j,it - emem1) = depe2d(i, j, nn2(it)) svor2(i,j,it - emem1) = epeflux2d(i, j, nn2(it)) stmp2(i,j,it - emem1) = baroc2d(i, j, nn2(it)) sspfh2(i,j,it - emem1) = diab2d(i, j, nn2(it)) shgt2(i,j,it - emem1) = corre2d(i, j, nn2(it)) saetoke2(i,j,it - emem1) = aetoke2d(i, j, nn2(it)) endif enddo enddo enddo !do it=1,itmax !write(6,*) sreh(20,20,nn2(it)),ssreh2(20,20,it),it !enddo do it=1,emem1 do j=jmn,jmx do i=imn,imx aveu1(i,j)=aveu1(i,j)+su1(i,j,it)/real(emem1) avev1(i,j)=avev1(i,j)+sv1(i,j,it)/real(emem1) avevor1(i,j)=avevor1(i,j)+svor1(i,j,it)/real(emem1) avetmp1(i,j)=avetmp1(i,j)+stmp1(i,j,it)/real(emem1) avespfh1(i,j)=avespfh1(i,j)+sspfh1(i,j,it)/real(emem1) avehgt1(i,j)=avehgt1(i,j)+shgt1(i,j,it)/real(emem1) aveaetoke1(i,j)=aveaetoke1(i,j)+saetoke1(i,j,it)/real(emem1) enddo enddo enddo do it=1,emem-emem1+1 do j=jmn,jmx do i=imn,imx aveu2(i,j)=aveu2(i,j)+su2(i,j,it)/real(emem - emem1 + 1) avev2(i,j)=avev2(i,j)+sv2(i,j,it)/real(emem - emem1 + 1) avevor2(i,j)=avevor2(i,j)+svor2(i,j,it)/real(emem - emem1 + 1) avetmp2(i,j)=avetmp2(i,j)+stmp2(i,j,it)/real(emem - emem1 + 1) avespfh2(i,j)=avespfh2(i,j)+sspfh2(i,j,it)/real(emem - emem1 + 1) avehgt2(i,j)=avehgt2(i,j)+shgt2(i,j,it)/real(emem - emem1 + 1) aveaetoke2(i,j)=aveaetoke2(i,j)+saetoke2(i,j,it)/real(emem - emem1 + 1) enddo enddo enddo do j=jmn,jmx do i=imn,imx du(i,j) = aveu1(i,j) - aveu2(i,j) dv(i,j) = avev1(i,j) - avev2(i,j) dvor(i,j) = avevor1(i,j) - avevor2(i,j) dtmp(i,j) = avetmp1(i,j) - avetmp2(i,j) dspfh(i,j) = avespfh1(i,j) - avespfh2(i,j) dhgt(i,j) = avehgt1(i,j) - avehgt2(i,j) daetoke(i,j) = aveaetoke1(i,j) - aveaetoke2(i,j) if(abs(du(i,j)).ge.abs(dcou(i,j))) then nu(i,j) = nu(i,j) + 1 endif if(abs(dv(i,j)).ge.abs(dcov(i,j))) then nv(i,j) = nv(i,j) + 1 endif if(abs(dvor(i,j)).ge.abs(dcovor(i,j))) then nvor(i,j) = nvor(i,j) + 1 endif if(abs(dtmp(i,j)).ge.abs(dcotmp(i,j))) then ntmp(i,j) = ntmp(i,j) + 1 endif if(abs(dspfh(i,j)).ge.abs(dcospfh(i,j))) then nspfh(i,j) = nspfh(i,j) + 1 endif if(abs(dhgt(i,j)).ge.abs(dcohgt(i,j))) then nhgt(i,j) = nhgt(i,j) + 1 endif if(abs(daetoke(i,j)).ge.abs(dcoaetoke(i,j))) then naetoke(i,j) = naetoke(i,j) + 1 endif enddo enddo !write(6,*) avesreh1(20,20),avesreh2(20,20),dsreh(20,20),dcosreh(20,20) !write(6,*) nnt enddo !nntloop do j=jmn,jmx do i=imn,imx pu(i,j) = (1.+real(nu(i,j)))/(real(nmx+1)) pv(i,j) = (1.+real(nv(i,j)))/(real(nmx+1)) pvor(i,j) = (1.+real(nvor(i,j)))/(real(nmx+1)) ptmp(i,j) = (1.+real(ntmp(i,j)))/(real(nmx+1)) pspfh(i,j) = (1.+real(nspfh(i,j)))/(real(nmx+1)) phgt(i,j) = (1.+real(nhgt(i,j)))/(real(nmx+1)) paetoke(i,j) = (1.+real(naetoke(i,j)))/(real(nmx+1)) enddo enddo !psreh=psreh*2. !pcold=pcold*2. !pwarm=pwarm*2. write(6,*) pu(20,20),pv(10,10),pvor(20,20) write(69) ((pu(i,j), i = imn,imx), j = jmn,jmx) write(69) ((pv(i,j), i = imn,imx), j = jmn,jmx) write(69) ((pvor(i,j),i = imn,imx), j = jmn,jmx) write(69) ((ptmp(i,j),i = imn,imx), j = jmn,jmx) write(69) ((pspfh(i,j),i = imn,imx), j = jmn,jmx) write(69) ((phgt(i,j),i = imn,imx), j = jmn,jmx) write(69) ((paetoke(i,j),i = imn,imx), j = jmn,jmx) endprogram permutation