da_print_be_stats_p.inc
References to this file elsewhere.
1 subroutine da_print_be_stats_p(outunit, ni, nj, nk, num_bins, num_bins2d, bin, &
2 bin2d, regcoeff1, regcoeff2, regcoeff3)
3
4 !----------------------------------------------------------------------------
5 ! Purpose: To print out the regression coefficients for the physical
6 ! transform.
7 !
8 ! Input : outunit --- Fortran unit for writing out
9 ! ni,nj,nk --- Dimensions
10 ! num_bins, num_bins2d --- Number of the 3d and 2d bins
11 ! bin, bin2d --- bin index for the gridpoints
12 ! regcoeff1 --- Reg coefs for chi = regcoeff1 * psi
13 ! regcoeff2 --- Reg coefs for Ps = sum(regcoeff2(1:k)*psi(1:k))
14 ! regcoeff3 --- Reg coefs for T(k)= sum(regcoeff3(k,1:k)*psi(1:k))
15 !
16 ! Output : fort.171 --- regcoef1
17 ! fort.172 --- regcoeff2
18 ! fort.173 --- regcoeff3
19 !----------------------------------------------------------------------------
20
21 implicit none
22
23 integer, intent(inout) :: outunit ! Output file unit.
24 integer, intent(in) :: ni, nj, nk ! Grid dimensions.
25 integer, intent(in) :: num_bins ! Number of 3D field bins.
26 integer, intent(in) :: num_bins2d ! Number of 2D field bins.
27 integer, intent(in) :: bin(1:ni,1:nj,1:nk) ! Bin assigned to each 3D point.
28 integer, intent(in) :: bin2d(1:ni,1:nj) ! Bin assigned to each 3D point.
29 real, intent(in) :: regcoeff1(1:num_bins) ! psi/chi regression cooefficient.
30 real, intent(in) :: regcoeff2(1:nk,1:num_bins2d) ! psi/ps regression cooefficient.
31 real, intent(in) :: regcoeff3(1:nk,1:nk,1:num_bins2d) ! psi/T regression cooefficient.
32
33 integer :: k1, k2, i, k, j, b, number ! Loop counters.
34
35 write(unit=stdout,fmt='(a,i5)')' psi/chi regression coefficient in unit ', outunit
36
37 open(unit=outunit)
38 write(unit=outunit,fmt='(a,i5)')' psi/chi regression coefficient in unit ', outunit
39 do b = 1, num_bins
40 number = 0
41 do i = 1,ni
42 do j = 1,nj
43 do k = 1,nk
44 if (bin(i,j,k) == b) then
45 number = number + 1
46 end if
47 end do
48 end do
49 end do
50 write(unit=outunit,fmt='("bin=",i6," number_in_bin=",i6,5X,1pE12.5)') &
51 b, number, regcoeff1(b)
52 end do
53 close(unit=outunit)
54
55 outunit = outunit + 1
56
57 write(unit=stdout,fmt='(a,i5)') &
58 ' psi/ps regression coefficient in unit ', outunit
59
60 open(unit=outunit)
61 write(unit=outunit,fmt='(a,i5)') &
62 ' psi/ps regression coefficient in unit ', outunit
63
64 do b = 1, num_bins2d
65 number = 0
66 do i = 1,ni
67 do j = 1,nj
68 if (bin2d(i,j) == b) then
69 number = number + 1
70 end if
71 end do
72 end do
73 write(unit=outunit,fmt='(/"bin=",i6," number_in_bin=",i6)') b, number
74 write(unit=outunit,fmt='((2X,9(i3,2x,1pE12.5)))') (k,regcoeff2(k,b), k=1,nk)
75 end do
76 close(unit=outunit)
77
78 outunit = outunit + 1
79
80 write(unit=stdout,fmt='(a,i5)') &
81 ' psi/T regression coefficient in unit ', outunit
82
83 open(unit=outunit)
84 write(unit=outunit,fmt='(a,i5)') &
85 ' psi/T regression coefficient in unit ', outunit
86 do b = 1, num_bins2d
87 number = 0
88 do i = 1,ni
89 do j = 1,nj
90 if (bin2d(i,j) == b) then
91 number = number + 1
92 end if
93 end do
94 end do
95
96 write(unit=outunit,fmt='(/"bin=",i6," number_in_bin=",i6)') b, number
97 do k1 = 1,nk
98 write(unit=outunit,fmt='(2X,"Temperature at k1=",i3)') k1
99 write(unit=outunit,fmt='((2X,9(i3,2x,1pE12.5)))') (k2,regcoeff3(k1,k2,b), k2=1,nk)
100 end do
101 end do
102 close(unit=outunit)
103
104 write(unit=stdout,fmt=*) ' '
105
106 outunit = outunit + 1
107
108 end subroutine da_print_be_stats_p
109
110