program falsedetection implicit none integer::ii,jj,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,tmx1,tmx2,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::day1,hour1,hmx,hmn,hhmn,hhmx,hmx2,ob,nob integer::k,ememax,emem,emem1 parameter(imx=51,jmx=51,hmx=37,syear=1995,eyear=2012,smonth=3,emonth=5,hmn=1,jmn=1,imn=1,tmn=1,hhmn=1,hhmx=40,hmx2=27,nntmx=10000,tmx1=2172,tmx2=1909,itmax=5000,nmx=10000,ememax=2000) integer,dimension(itmax)::nn2 real,dimension(itmax)::r,rr real,dimension(nntmx)::trid3,time3 real,dimension(jmx)::co,lat real,dimension(jmx,itmax)::sco3,sco4,sco5 real,dimension(imx,jmx)::psreh, pbaroc, pdiab real,dimension(imx,jmx)::pu, pv, pvor, ptmp, pspfh, ppv, phgt real,dimension(imx,jmx,hmx2,itmax)::ssp1,ssp2 real::snum,snum1,snum2,snum3,tnum,tnum2 real, dimension(1000) :: pfdr, minsreh3, minbaroc3, mindiab3 integer, dimension(1000) :: isreh3, jsreh3 integer, dimension(1000) :: ibaroc3, jbaroc3 integer, dimension(1000) :: idiab3, jdiab3 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,jetpos*8 character:: jetpos1*8, jetpos2*8 logical::fast integer::imn2,jmn2,imx2,jmx2 integer::nmax,isreh,jsreh,ndata real :: alpha, minpsreh, minpsreh2, lonsreh, latsreh real :: minpbaroc, minpbaroc2, lonbaroc, latbaroc real :: minpdiab, minpdiab2, londiab, latdiab CAT='PO' CAT2='PO' jetpos1='Njetconf' jetpos2='Njetdiff' SEA='DJF' region='JPN' it=7 do it = 9, 9, 2 if(it.le.9) then write(kt1,'(i1)') it kt='0'//kt1 else write(kt,'(i2)') it endif !open(69,file='pvalue_for_permutation_frontxy_KT'//kt//'_'//CAT//'_'//CAT2//'_'//jetpos//'.data',status='unknown',form='unformatted') open(69,file='pvalue_for_permutation_ecjet_fgeneenv_rotate_'//CAT//'_'//kt//'_'//jetpos1//'_'//jetpos2//'.data',status='old',form='unformatted') open(70,file='false_detection_for_permutation_fgeneenv_def_rotate_850hPa_KT'//kt//'_'//CAT//'_'//jetpos1//'_'//jetpos2//'.txt',status='unknown',form='formatted') open(71,file='false_detection_for_permutation_fgeneenv_slant_rotate_850hPa_KT'//kt//'_'//CAT//'_'//jetpos1//'_'//jetpos2//'.txt',status='unknown',form='formatted') open(72,file='false_detection_for_permutation_fgeneenv_pt_rotate_850hPa_KT'//kt//'_'//CAT//'_'//jetpos1//'_'//jetpos2//'.txt',status='unknown',form='formatted') open(73,file='false_detection_for_permutation_fgeneenv_fmag_rotate_850hPa_KT'//kt//'_'//CAT//'_'//jetpos1//'_'//jetpos2//'.txt',status='unknown',form='formatted') ! read data read(69) ((pu(i, j), i = imn, imx), j = jmx, jmn, -1) read(69) ((pv(i, j), i = imn, imx), j = jmx, jmn, -1) read(69) ((pvor(i, j), i = imn, imx), j = jmx, jmn, -1) read(69) ((ptmp(i,j), i = imn, imx),j = jmx, jmn, -1) !write(6,*) psreh ! sort pvalues !imn2=10 !imx2=42 !jmn2=10 !jmx2=42 !imn2=9 !imx2=43 !jmn2=9 !jmx2=43 imn2=18 imx2=34 jmn2=18 jmx2=34 !imn2=14 !imx2=38 !jmn2=14 !jmx2=38 call fdr(imx, jmx, imn2, imx2, jmn2, jmx2, pu, ibaroc3, jbaroc3, nmax, pfdr, minbaroc3) write(6, *) nmax do n = 1, nmax lonbaroc = 230.+ (ibaroc3(n)-1) * 1.5 latbaroc = 10.+ (jbaroc3(n)-1) * 1.5 !write(70,*) minsreh3(n),isreh3(n),jsreh3(n),pfdr(n) write(70,*) lonbaroc, latbaroc, pfdr(n) write(6,*) minbaroc3(n), lonbaroc, latbaroc, pfdr(n) enddo nmax = 0 pfdr = 0. call fdr(imx, jmx, imn2, imx2, jmn2, jmx2, pv, idiab3, jdiab3, nmax, pfdr, mindiab3) write(6, *) nmax do n = 1, nmax londiab = 230. + (idiab3(n)-1) * 1.5 latdiab = 10. + (jdiab3(n)-1) * 1.5 !write(70,*) minsreh3(n),isreh3(n),jsreh3(n),pfdr(n) write(71,*) londiab, latdiab, pfdr(n) write(6,*) mindiab3(n), londiab, latdiab, pfdr(n) enddo nmax = 0 pfdr = 0. call fdr(imx, jmx, imn2, imx2, jmn2, jmx2, pvor, idiab3, jdiab3, nmax, pfdr, mindiab3) write(6, *) nmax do n = 1, nmax londiab = 230. + (idiab3(n)-1) * 1.5 latdiab = 10. + (jdiab3(n)-1) * 1.5 !write(70,*) minsreh3(n),isreh3(n),jsreh3(n),pfdr(n) write(72,*) londiab, latdiab, pfdr(n) write(6,*) mindiab3(n), londiab, latdiab, pfdr(n) enddo nmax = 0 pfdr = 0. call fdr(imx, jmx, imn2, imx2, jmn2, jmx2, ptmp, idiab3, jdiab3, nmax, pfdr, mindiab3) write(6, *) nmax do n = 1, nmax londiab = 230. + (idiab3(n)-1) * 1.5 latdiab = 10. + (jdiab3(n)-1) * 1.5 !write(70,*) minsreh3(n),isreh3(n),jsreh3(n),pfdr(n) write(73,*) londiab, latdiab, pfdr(n) write(6,*) mindiab3(n), londiab, latdiab, pfdr(n) enddo write(6,*) it write(6,*) (jmx2-jmn2+1)*(imx2-imn2+1) enddo endprogram falsedetection subroutine fdr(imx, jmx, imn2, imx2, jmn2, jmx2, var, ivar3, jvar3, nmax, pfdr, minvar3) integer :: i, j, ii, jj integer :: nmax integer, intent(in) :: imx, jmx integer, intent(in) :: imn2, imx2, jmn2, jmx2 real, intent(in), dimension(imx, jmx) :: var integer :: n integer :: num, ndata real :: alpha real :: minvar, minvar2 integer :: ivar, jvar integer, intent(out), dimension(1000) :: ivar3, jvar3 real, intent(out), dimension(1000) :: minvar3, pfdr ndata=(jmx2 - jmn2 + 1)*(imx2 - imn2 + 1) !write(6, *) var ! sort pvalues alpha = 0.1 num=0 minvar2 = -1000 do jj = jmn2, jmx2 do ii = imn2, imx2 minvar = 1. do j = jmn2, jmx2 do i = imn2, imx2 if(var(i, j) .gt. minvar2 .and. var(i, j) .lt. minvar) then minvar = var(i,j) ivar = i jvar = j endif enddo enddo num = num+1 minvar2 = minvar minvar3(num) = minvar ivar3(num) = ivar jvar3(num) = jvar pfdr(num) = num*alpha/real(ndata) ! write(6,*) minvar3(num),ivar3(num),jvar3(num),pfdr(num) enddo enddo nmax = 0 do n = 1, num if( minvar3(n) .lt. pfdr(n) .and. n.gt.nmax ) then nmax = n endif enddo write(6, *) num end subroutine fdr