source: codes/icosagcm/devel/src/vertical/vertical_remap.f90

Last change on this file was 961, checked in by jisesh, 5 years ago

devel: fixed array index in vertical_remap

File size: 4.9 KB
Line 
1MODULE vertical_remap_mod
2  USE icosa
3  USE omp_para
4  IMPLICIT NONE
5  PRIVATE
6PUBLIC vertical_remap_extdata,compute_vertical_remap_extdata
7
8  PUBLIC :: vertical_remap
9
10CONTAINS
11
12   
13  SUBROUTINE vertical_remap(pressure_level,field_in,f_ps,field_out)
14    USE compute_pression_mod, ONLY : pression
15    REAL(rstd), INTENT(IN) :: pressure_level(:)
16    TYPE(t_field),POINTER :: field_in(:)
17    TYPE(t_field),POINTER :: f_ps(:)
18    TYPE(t_field),POINTER :: field_out(:)
19
20    TYPE(t_field),POINTER,SAVE :: f_p(:)
21    REAL(rstd),POINTER :: in(:,:)
22    REAL(rstd),POINTER :: out(:,:)
23    REAL(rstd),POINTER :: p(:,:)
24   
25    INTEGER :: ind
26
27    CALL allocate_field(f_p,field_t,type_real,llm+1)
28    CALL pression(f_ps,f_p)
29 
30    DO ind=1,ndomain
31      IF (.NOT. assigned_domain(ind)) CYCLE
32      CALL swap_dimensions(ind)
33      CALL swap_geometry(ind)
34      p=f_p(ind)
35      in=field_in(ind)
36      out=field_out(ind)
37      CALL compute_vertical_remap(pressure_level,in,p,out)
38    ENDDO
39   
40  END SUBROUTINE  vertical_remap
41
42  SUBROUTINE compute_vertical_remap(pressure_level,in,p,out)
43    REAL(rstd),INTENT(IN)  :: pressure_level(:)
44    REAL(rstd),INTENT(IN)  :: in(:,:)
45    REAL(rstd),INTENT(IN)  :: p(iim*jjm,llm+1)
46    REAL(rstd),INTENT(OUT) :: out(iim*jjm,llm)
47    REAL(rstd) :: coeff, pmid
48    INTEGER :: i,j,ij,l,n,nb_level
49    INTEGER :: a
50    INTEGER :: b
51    LOGICAL :: positive
52   
53    nb_level=size(pressure_level)
54    IF (pressure_level(1)>=pressure_level(nb_level)) THEN
55      positive=.FALSE.
56    ELSE
57      positive=.TRUE.
58    ENDIF
59     
60 !$OMP BARRIER   
61    IF (is_omp_level_master) THEN
62   
63    DO l=1,llm
64      DO j=jj_begin,jj_end
65        DO i=ii_begin,ii_end
66          ij=(j-1)*iim+i
67          pmid=0.5*(p(ij,l)+p(ij,l+1))
68          IF (positive) THEN
69            a=0
70            DO n=1,nb_level-1
71              IF ( (pmid>=pressure_level(n) .AND. pmid<pressure_level(n+1))) THEN
72               a=n ; b=n+1 ; EXIT
73              ENDIF
74            ENDDO
75            IF (a==0) THEN
76              IF (pmid<=pressure_level(1)) THEN
77                a=1 ; b=2
78              ELSE
79                a=nb_level-1 ; b=nb_level
80              ENDIF
81            ENDIF
82          ELSE
83            a=0
84            DO n=1,nb_level-1
85              IF ( (pmid<=pressure_level(n) .AND. pmid>pressure_level(n+1))) THEN
86               a=n ; b=n+1 ; EXIT
87              ENDIF
88            ENDDO
89           
90            IF (a==0) THEN
91              IF (pmid >= pressure_level(1)) THEN
92                a=1 ; b=2
93              ELSE
94                a=nb_level-1 ; b=nb_level
95              ENDIF
96            ENDIF
97          ENDIF
98                 
99          coeff=(pmid-pressure_level(a))/(pressure_level(a)-pressure_level(b))
100          out(ij,l)=in(ij,a)+coeff*(in(ij,a)-in(ij,b))
101        ENDDO
102      ENDDO
103    ENDDO
104 
105    ENDIF
106 !$OMP BARRIER   
107       
108  END SUBROUTINE compute_vertical_remap
109
110  SUBROUTINE vertical_remap_extdata(field_in,f_target_pressure,field_out)
111  USE icosa
112  USE omp_para
113  USE disvert_mod, ONLY : presnivs
114
115  IMPLICIT NONE
116    TYPE(t_field),POINTER :: field_in(:) 
117    TYPE(t_field),POINTER :: field_out(:)
118    TYPE(t_field),POINTER :: f_target_pressure(:)
119
120    REAL(rstd),POINTER :: target_pressure(:,:)
121    REAL(rstd),POINTER :: in(:,:)
122    REAL(rstd),POINTER :: out(:,:)
123    INTEGER :: ind
124
125    DO ind=1,ndomain
126      IF (.NOT. assigned_domain(ind)) CYCLE
127      CALL swap_dimensions(ind)
128      CALL swap_geometry(ind)
129      in=field_in(ind)
130      out=field_out(ind)
131      target_pressure=f_target_pressure(ind)
132      CALL compute_vertical_remap_extdata(in,target_pressure,out)
133    ENDDO
134
135  END SUBROUTINE vertical_remap_extdata
136
137  SUBROUTINE compute_vertical_remap_extdata(in,target_pressure,out2d_press)
138  USE omp_para
139  USE disvert_mod, ONLY : presnivs 
140  IMPLICIT NONE
141    REAL(rstd),INTENT(IN)  :: in(:,:)
142    REAL(rstd),INTENT(IN)  :: target_pressure(iim*jjm,llm) 
143    REAL(rstd),INTENT(OUT) :: out2d_press(iim*jjm,llm)
144    REAL(rstd) :: coeff, target_pval,testp1,testp2
145    INTEGER :: i,j,ij,l,n,nb_level
146    INTEGER :: a
147    INTEGER :: b
148    LOGICAL :: positive
149   
150    nb_level=size(presnivs)
151 !$OMP BARRIER   
152    IF (is_omp_level_master) THEN
153    DO l=1,llm
154      DO j=jj_begin,jj_end
155        DO i=ii_begin,ii_end
156          ij=(j-1)*iim+i
157          target_pval=target_pressure(ij,l)
158            a=0
159            DO n=1,nb_level-1
160              IF ( (target_pval<=presnivs(n) .AND. target_pval>=presnivs(n+1))) THEN
161               testp1=presnivs(n); testp2=presnivs(n+1)
162               a=n ; b=n+1 ; EXIT
163              ENDIF
164            ENDDO
165            IF (a==0) THEN
166              IF (target_pval>=presnivs(1)) THEN
167                a=1 ; b=2
168              ELSE
169                a=nb_level-1 ; b=nb_level
170              ENDIF
171            ENDIF
172
173          coeff=(target_pval-presnivs(a))/(presnivs(a)-presnivs(b))
174          out2d_press(ij,l)=in(ij,a)+coeff*(in(ij,a)-in(ij,b))
175        ENDDO
176      ENDDO
177    ENDDO
178
179    ENDIF
180 !$OMP BARRIER   
181
182  END SUBROUTINE compute_vertical_remap_extdata
183
184
185
186END MODULE vertical_remap_mod
Note: See TracBrowser for help on using the repository browser.