Use NetCDF
     Real, Allocatable :: A(:,:), B(:,:), mode(:), time(:), X(:), Y(:)
     Real, Parameter :: Miss=-99.0
     Integer :: i,j,k, Nmode, Ntime
     Real :: Corr,prob
     Character (Len=255) :: Fname

     Fname="/home/tanpv/Progs/2021/NC/AB_Dr_SST.nc"
     Call rd_ncdf_2d (Fname,mode,time,Nmode,Ntime,A,"a")
     Deallocate (mode,time)
     Call rd_ncdf_2d (Fname,mode,time,Nmode,Ntime,B,"b")
     Allocate (X(Ntime),Y(Ntime))
     Do i=1,Nmode
        X(:)=A(i,:); Y(:)=B(i,:)
        Call Correlation(X,Y,Ntime,Corr,prob,Miss)
        print*,i,Ntime,Corr,prob
     Enddo

     CONTAINS

      !--------------------------------------------------------------------
      Subroutine Correlation(X,Y,N,Corr,prob,Miss)
        Integer :: N, N1, i, j
        Real :: X(n), Y(n), Miss, Xtb, Ytb, Sx, Sy, Rxy, tmp, t

        N1=0; Xtb=0; Ytb=0; Sx=0; Sy=0; Rxy=0
        Do i=1,N
           if (X(i)/=Miss.AND.Y(i)/=Miss) then
              n1=n1+1
              Xtb=Xtb+X(i); Sx=Sx+X(i)*X(i)
              Ytb=Ytb+Y(i); Sy=Sy+Y(i)*Y(i)
              Rxy=Rxy+X(i)*Y(i)
           endif
        Enddo
        if (N1>2) then
           Xtb=Xtb/Real(N1); Ytb=Ytb/Real(N1)
           Sx=Sx/Real(N1)-Xtb*Xtb; Sy=Sy/Real(N1)-Ytb*Ytb
           tmp=Rxy/Real(N1)
           Rxy=(tmp-Xtb*Ytb)/Sqrt(Sx*Sy)
           t = ABS(Rxy/ ( Sqrt(1-Rxy**2) / Sqrt(Real(N1-2)) ) )
           Prob = T_Dist (t, Real(N1-2))
        else
           Rxy=Miss; prob=Miss
        endif
        Corr=Rxy

      End Subroutine Correlation
      !---------------------------------------------------------------------------
      Include "/home/tanpv/Mylibs/sub_rd_ncdf_Tan.F90"

      END
      Include "/home/tanpv/Mylibs/Student_Prob.f90"