Context Navigation

module_ForDiagnostics.f90 @ 2210

Last change on this file since 2210 was 2209, checked in by lfita, 7 years ago
Adding `range_faces'
File size: 34.4 KB

Rev	Line
[770]	1	!! Fortran version of different diagnostics
	2	! L. Fita. LMD May 2016
[772]	3	! gfortran module_generic.o module_ForDiagnosticsVars.o -c module_ForDiagnostics.F90
	4	!
	5	! f2py -m module_ForDiagnostics --f90exec=/usr/bin/gfortran-4.7 -c module_generic.F90 module_ForDiagnosticsVars.F90 module_ForDiagnostics.F90
	6
[770]	7	MODULE module_ForDiagnostics
	8
[1608]	9	USE module_definitions
[770]	10	USE module_generic
[772]	11	USE module_ForDiagnosticsVars
[770]	12
[772]	13	CONTAINS
[770]	14
[772]	15	!!!!!!! Calculations
[1769]	16	! compute_cape_afwa4D: Subroutine to use WRF phys/module_diag_afwa.F `buyoancy' subroutine to compute
	17	! CAPE, CIN, ZLFC, PLFC, LI
[1758]	18	! compute_cllmh4D3: Computation of low, medium and high cloudiness from a 4D CLDFRA and pressure being
	19	! 3rd dimension the z-dim
	20	! compute_cllmh3D3: Computation of low, medium and high cloudiness from a 3D CLDFRA and pressure being
	21	! 3rd dimension the z-dim
[772]	22	! compute_cllmh: Computation of low, medium and high cloudiness
	23	! compute_clt4D3: Computation of total cloudiness from a 4D CLDFRA being 3rd dimension the z-dim
	24	! compute_clt3D3: Computation of total cloudiness from a 3D CLDFRA being 3rd dimension the z-dim
	25	! compute_clt: Computation of total cloudiness
[1908]	26	! compute_fog_K84: Computation of fog and visibility following Kunkel, (1984)
	27	! compute_fog_RUC: Computation of fog and visibility following RUC method Smirnova, (2000)
[1909]	28	! compute_fog_FRAML50: fog and visibility following Gultepe and Milbrandt, (2010)
[2209]	29	! compute_massvertint1D: Subroutine to vertically integrate a 1D variable in eta vertical coordinates
[1795]	30	! compute_psl_ecmwf: Compute sea level pressure using ECMWF method following Mats Hamrud and Philippe Courtier [Pa]
[2209]	31	! compute_range_faces: Subroutine to compute faces [uphill, valleys, downhill] of a mountain range along a given face
[1769]	32	! compute_vertint1D: Subroutine to vertically integrate a 1D variable in any vertical coordinates
	33	! compute_zint4D: Subroutine to vertically integrate a 4D variable in any vertical coordinates
[1773]	34	! compute_zmla_generic4D: Subroutine to compute pbl-height following a generic method
[1776]	35	! compute_zwind4D: Subroutine to compute extrapolate the wind at a given height following the 'power law' methodology
[1784]	36	! compute_zwind_log4D: Subroutine to compute extrapolate the wind at a given height following the 'logarithmic law' methodology
[1783]	37	! compute_zwindMCO3D: Subroutine to compute extrapolate the wind at a given height following the 'power law' methodolog
[1773]	38
[772]	39	!!!
	40	! Calculations
	41	!!!
[770]	42
[772]	43	SUBROUTINE compute_cllmh4D2(cldfra4D, pres4D, cllmh4D2, d1, d2, d3, d4)
	44	! Subroutine to compute the low, medium and high cloudiness following 'newmicro.F90' from LMDZ from a 4D CLDFRA and pressure
	45	! where zdim is the 2nd dimension (thus, cldfra4D(d1,d2,d3,d4) --> cllmh(3,d1,d3,d4) 1: low, 2: medium, 3: high
	46	! It should be properly done via an 'INTERFACE', but...
[770]	47
[772]	48	IMPLICIT NONE
	49
[1141]	50	INTEGER, INTENT(in) :: d1, d2, d3, d4
[772]	51	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: cldfra4D, pres4D
	52	REAL(r_k), DIMENSION(3,d1,d3,d4), INTENT(out) :: cllmh4D2
	53
	54	! Local
	55	INTEGER :: i,j,k, zdim, Ndim
	56
[770]	57	!!!!!!! Variables
[772]	58	! cldfra4D: 4D cloud fraction values [1]
	59	! pres4D: 4D pressure values [Pa]
	60	! Ndim: number of dimensions of the input data
	61	! d[1-4]: dimensions of 'cldfra4D'
	62	! zdim: number of the vertical-dimension within the matrix
	63	! cltlmh4D2: low, medium, high cloudiness for the 4D cldfra and d2 being 'zdim'
[770]	64
[772]	65	fname = 'compute_cllmh4D2'
	66	zdim = 2
	67	Ndim = 4
[770]	68
[772]	69	DO i=1, d1
	70	DO j=1, d3
	71	DO k=1, d4
	72	cllmh4D2(:,i,j,k) = var_cllmh(cldfra4D(i,:,j,k), pres4D(i,:,j,k), d2)
	73	END DO
	74	END DO
	75	END DO
	76
	77	RETURN
[770]	78
[772]	79	END SUBROUTINE compute_cllmh4D2
[770]	80
[772]	81	SUBROUTINE compute_cllmh3D1(cldfra3D, pres3D, cllmh3D1, d1, d2, d3)
	82	! Subroutine to compute the low, medium and high cloudiness following 'newmicro.F90' from LMDZ from a 3D CLDFRA and pressure
	83	! where zdim is the 1st dimension (thus, cldfra3D(d1,d2,d3) --> cllmh(3,d2,d3) 1: low, 2: medium, 3: high
	84	! It should be properly done via an 'INTERFACE', but...
[770]	85
[772]	86	IMPLICIT NONE
[770]	87
[1141]	88	INTEGER, INTENT(in) :: d1, d2, d3
[772]	89	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: cldfra3D, pres3D
	90	REAL(r_k), DIMENSION(3,d2,d3), INTENT(out) :: cllmh3D1
	91
	92	! Local
	93	INTEGER :: i,j,k, zdim, Ndim
	94
	95	!!!!!!! Variables
	96	! cldfra3D: 3D cloud fraction values [1]
	97	! pres3D: 3D pressure values [Pa]
	98	! Ndim: number of dimensions of the input data
	99	! d[1-3]: dimensions of 'cldfra3D'
	100	! zdim: number of the vertical-dimension within the matrix
	101	! cltlmh3D1: low, medium, high cloudiness for the 3D cldfra and d1 being 'zdim'
	102
	103	fname = 'compute_cllmh3D1'
	104	zdim = 1
	105	Ndim = 3
	106
	107	DO i=1, d1
	108	DO j=1, d2
	109	cllmh3D1(:,i,j) = var_cllmh(cldfra3D(:,i,j), pres3D(:,i,j), d1)
	110	END DO
	111	END DO
	112
	113	RETURN
	114
	115	END SUBROUTINE compute_cllmh3D1
	116
	117	SUBROUTINE compute_cllmh(cldfra1D, cldfra2D, cldfra3D, cldfra4D, pres1D, pres2D, pres3D, pres4D, &
	118	Ndim, zdim, cllmh1D, cllmh2D1, cllmh2D2, cllmh3D1, cllmh3D2, cllmh3D3, cllmh4D1, cllmh4D2, &
	119	cllmh4D3, cllmh4D4, d1, d2, d3, d4)
	120	! Subroutine to compute the low, medium and high cloudiness following 'newmicro.F90' from LMDZ
	121
[770]	122	IMPLICIT NONE
	123
[1141]	124	INTEGER, INTENT(in) :: Ndim, d1, d2, d3, d4, zdim
[772]	125	REAL(r_k), DIMENSION(d1), OPTIONAL, INTENT(in) :: cldfra1D, pres1D
	126	REAL(r_k), DIMENSION(d1,d2), OPTIONAL, INTENT(in) :: cldfra2D, pres2D
	127	REAL(r_k), DIMENSION(d1,d2,d3), OPTIONAL, INTENT(in) :: cldfra3D, pres3D
	128	REAL(r_k), DIMENSION(d1,d2,d3,d4), OPTIONAL, &
	129	INTENT(in) :: cldfra4D, pres4D
	130	REAL(r_k), DIMENSION(3), OPTIONAL, INTENT(out) :: cllmh1D
	131	REAL(r_k), DIMENSION(d1,3), OPTIONAL, INTENT(out) :: cllmh2D1
	132	REAL(r_k), DIMENSION(d2,3), OPTIONAL, INTENT(out) :: cllmh2D2
	133	REAL(r_k), DIMENSION(d2,d3,3), OPTIONAL, INTENT(out) :: cllmh3D1
	134	REAL(r_k), DIMENSION(d1,d3,3), OPTIONAL, INTENT(out) :: cllmh3D2
	135	REAL(r_k), DIMENSION(d1,d2,3), OPTIONAL, INTENT(out) :: cllmh3D3
	136	REAL(r_k), DIMENSION(d2,d3,d4,3), OPTIONAL, &
	137	INTENT(out) :: cllmh4D1
	138	REAL(r_k), DIMENSION(d1,d3,d4,3), OPTIONAL, &
	139	INTENT(out) :: cllmh4D2
	140	REAL(r_k), DIMENSION(d1,d2,d4,3), OPTIONAL, &
	141	INTENT(out) :: cllmh4D3
	142	REAL(r_k), DIMENSION(d1,d2,d3,3), OPTIONAL, &
	143	INTENT(out) :: cllmh4D4
[770]	144
	145	! Local
[772]	146	INTEGER :: i,j,k
[770]	147
	148	!!!!!!! Variables
[772]	149	! cldfra[1-4]D: cloud fraction values [1]
	150	! pres[1-4]D: pressure values [Pa]
	151	! Ndim: number of dimensions of the input data
	152	! d[1-4]: dimensions of 'cldfra'
	153	! zdim: number of the vertical-dimension within the matrix
	154	! cllmh1D: low, medium and high cloudiness for the 1D cldfra
	155	! cllmh2D1: low, medium and high cloudiness for the 2D cldfra and d1 being 'zdim'
	156	! cllmh2D2: low, medium and high cloudiness for the 2D cldfra and d2 being 'zdim'
	157	! cllmh3D1: low, medium and high cloudiness for the 3D cldfra and d1 being 'zdim'
	158	! cllmh3D2: low, medium and high cloudiness for the 3D cldfra and d2 being 'zdim'
	159	! cllmh3D3: low, medium and high cloudiness for the 3D cldfra and d3 being 'zdim'
	160	! cllmh4D1: low, medium and high cloudiness for the 4D cldfra and d1 being 'zdim'
	161	! cllmh4D2: low, medium and high cloudiness for the 4D cldfra and d2 being 'zdim'
	162	! cllmh4D3: low, medium and high cloudiness for the 4D cldfra and d3 being 'zdim'
	163	! cllmh4D4: low, medium and high cloudiness for the 4D cldfra and d4 being 'zdim'
[770]	164
[772]	165	fname = 'compute_cllmh'
[770]	166
[772]	167	SELECT CASE (Ndim)
	168	CASE (1)
	169	cllmh1D = var_cllmh(cldfra1D, pres1D, d1)
	170	CASE (2)
	171	IF (zdim == 1) THEN
	172	DO i=1, d2
	173	cllmh2D1(i,:) = var_cllmh(cldfra2D(:,i), pres2D(:,i), d1)
	174	END DO
	175	ELSE IF (zdim == 2) THEN
	176	DO i=1, d1
	177	cllmh2D2(i,:) = var_cllmh(cldfra2D(:,i), pres2D(i,:), d2)
	178	END DO
	179	ELSE
	180	PRINT *,TRIM(ErrWarnMsg('err'))
	181	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	182	PRINT *,' accepted values: 1,2'
	183	STOP
	184	END IF
	185	CASE (3)
	186	IF (zdim == 1) THEN
	187	DO i=1, d2
	188	DO j=1, d3
	189	cllmh3D1(i,j,:) = var_cllmh(cldfra3D(:,i,j), pres3D(:,i,j), d1)
	190	END DO
	191	END DO
	192	ELSE IF (zdim == 2) THEN
	193	DO i=1, d1
	194	DO j=1, d3
	195	cllmh3D2(i,j,:) = var_cllmh(cldfra3D(i,:,j), pres3D(i,:,j), d2)
	196	END DO
	197	END DO
	198	ELSE IF (zdim == 3) THEN
	199	DO i=1, d1
	200	DO j=1, d2
	201	cllmh3D3(i,j,:) = var_cllmh(cldfra3D(i,j,:), pres3D(i,j,:), d3)
	202	END DO
	203	END DO
	204	ELSE
	205	PRINT *,TRIM(ErrWarnMsg('err'))
	206	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	207	PRINT *,' accepted values: 1,2,3'
	208	STOP
	209	END IF
	210	CASE (4)
	211	IF (zdim == 1) THEN
	212	DO i=1, d2
	213	DO j=1, d3
	214	DO k=1, d4
	215	cllmh4D1(i,j,k,:) = var_cllmh(cldfra4D(:,i,j,k), pres4D(:,i,j,k), d1)
	216	END DO
	217	END DO
	218	END DO
	219	ELSE IF (zdim == 2) THEN
	220	DO i=1, d1
	221	DO j=1, d3
	222	DO k=1, d4
	223	cllmh4D2(i,j,k,:) = var_cllmh(cldfra4D(i,:,j,k), pres4D(i,:,j,k), d2)
	224	END DO
	225	END DO
	226	END DO
	227	ELSE IF (zdim == 3) THEN
	228	DO i=1, d2
	229	DO j=1, d3
	230	DO k=1, d4
	231	cllmh4D3(i,j,k,:) = var_cllmh(cldfra4D(i,j,:,k), pres4D(i,j,:,k), d3)
	232	END DO
	233	END DO
	234	END DO
	235	ELSE IF (zdim == 4) THEN
	236	DO i=1, d1
	237	DO j=1, d2
	238	DO k=1, d3
	239	cllmh4D4(i,j,k,:) = var_cllmh(cldfra4D(i,j,k,:), pres4D(i,j,k,:), d4)
	240	END DO
	241	END DO
	242	END DO
	243	ELSE
	244	PRINT *,TRIM(ErrWarnMsg('err'))
	245	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	246	PRINT *,' accepted values: 1,2,3,4'
	247	STOP
	248	END IF
	249	CASE DEFAULT
	250	PRINT *,TRIM(ErrWarnMsg('err'))
	251	PRINT *,' ' // TRIM(fname) // ': Ndim:', Ndim,' not ready !!'
	252	STOP
	253	END SELECT
[770]	254
	255	RETURN
	256
[772]	257	END SUBROUTINE compute_cllmh
[770]	258
[772]	259	SUBROUTINE compute_clt4D2(cldfra4D, clt4D2, d1, d2, d3, d4)
	260	! Subroutine to compute the total cloudiness following 'newmicro.F90' from LMDZ from a 4D CLDFRA
	261	! where zdim is the 2nd dimension (thus, cldfra4D(d1,d2,d3,d4) --> clt(d1,d3,d4)
	262	! It should be properly done via an 'INTERFACE', but...
[770]	263
	264	IMPLICIT NONE
	265
[1141]	266	INTEGER, INTENT(in) :: d1, d2, d3, d4
[772]	267	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: cldfra4D
	268	REAL(r_k), DIMENSION(d1,d3,d4), INTENT(out) :: clt4D2
	269
[770]	270	! Local
[772]	271	INTEGER :: i,j,k, zdim, Ndim
	272
[770]	273	!!!!!!! Variables
[772]	274	! cldfra4D: 4D cloud fraction values [1]
	275	! Ndim: number of dimensions of the input data
	276	! d[1-4]: dimensions of 'cldfra4D'
	277	! zdim: number of the vertical-dimension within the matrix
	278	! clt4D2: total cloudiness for the 4D cldfra and d2 being 'zdim'
[770]	279
[772]	280	fname = 'compute_clt4D2'
	281	zdim = 2
	282	Ndim = 4
[770]	283
[772]	284	DO i=1, d1
	285	DO j=1, d3
	286	DO k=1, d4
	287	clt4D2(i,j,k) = var_clt(cldfra4D(i,:,j,k), d2)
	288	END DO
	289	END DO
[770]	290	END DO
[772]	291
	292	RETURN
[770]	293
[772]	294	END SUBROUTINE compute_clt4D2
[770]	295
[772]	296	SUBROUTINE compute_clt3D1(cldfra3D, clt3D1, d1, d2, d3)
	297	! Subroutine to compute the total cloudiness following 'newmicro.F90' from LMDZ from a 3D CLDFRA
	298	! where zdim is the 1st dimension (thus, cldfra4D(d1,d2,d3) --> clt(d2,d3)
	299	! It should be properly done via an 'INTERFACE', but...
[770]	300
[772]	301	IMPLICIT NONE
[770]	302
[1141]	303	INTEGER, INTENT(in) :: d1, d2, d3
[772]	304	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: cldfra3D
	305	REAL(r_k), DIMENSION(d2,d3), INTENT(out) :: clt3D1
[770]	306
[772]	307	! Local
	308	INTEGER :: i,j,k, zdim, Ndim
	309
	310	!!!!!!! Variables
	311	! cldfra3D: 3D cloud fraction values [1]
	312	! Ndim: number of dimensions of the input data
	313	! d[1-3]: dimensions of 'cldfra3D'
	314	! zdim: number of the vertical-dimension within the matrix
	315	! clt3D1: total cloudiness for the 3D cldfra and d1 being 'zdim'
	316
	317	fname = 'compute_clt3D1'
	318	zdim = 1
	319	Ndim = 3
	320
	321	DO i=1, d2
	322	DO j=1, d3
	323	clt3D1(i,j) = var_clt(cldfra3D(:,i,j), d1)
	324	END DO
	325	END DO
	326
	327	RETURN
	328
	329	END SUBROUTINE compute_clt3D1
	330
	331	SUBROUTINE compute_clt(cldfra1D, cldfra2D, cldfra3D, cldfra4D, Ndim, zdim, clt1D, clt2D1, clt2D2, &
	332	clt3D1, clt3D2, clt3D3, clt4D1, clt4D2, clt4D3, clt4D4, d1, d2, d3, d4)
	333	! Subroutine to compute the total cloudiness following 'newmicro.F90' from LMDZ
	334
[770]	335	IMPLICIT NONE
	336
[1141]	337	INTEGER, INTENT(in) :: Ndim, d1, d2, d3, d4, zdim
[770]	338	REAL(r_k), DIMENSION(d1), OPTIONAL, INTENT(in) :: cldfra1D
	339	REAL(r_k), DIMENSION(d1,d2), OPTIONAL, INTENT(in) :: cldfra2D
	340	REAL(r_k), DIMENSION(d1,d2,d3), OPTIONAL, INTENT(in) :: cldfra3D
	341	REAL(r_k), DIMENSION(d1,d2,d3,d4), OPTIONAL, &
	342	INTENT(in) :: cldfra4D
	343	REAL(r_k), OPTIONAL, INTENT(out) :: clt1D
	344	REAL(r_k), DIMENSION(d1), OPTIONAL, INTENT(out) :: clt2D1
	345	REAL(r_k), DIMENSION(d2), OPTIONAL, INTENT(out) :: clt2D2
	346	REAL(r_k), DIMENSION(d2,d3), OPTIONAL, INTENT(out) :: clt3D1
	347	REAL(r_k), DIMENSION(d1,d3), OPTIONAL, INTENT(out) :: clt3D2
	348	REAL(r_k), DIMENSION(d1,d2), OPTIONAL, INTENT(out) :: clt3D3
	349	REAL(r_k), DIMENSION(d2,d3,d4), OPTIONAL,INTENT(out) :: clt4D1
	350	REAL(r_k), DIMENSION(d1,d3,d4), OPTIONAL,INTENT(out) :: clt4D2
	351	REAL(r_k), DIMENSION(d1,d2,d4), OPTIONAL,INTENT(out) :: clt4D3
	352	REAL(r_k), DIMENSION(d1,d2,d3), OPTIONAL,INTENT(out) :: clt4D4
	353
	354	! Local
	355	INTEGER :: i,j,k
	356
	357	!!!!!!! Variables
	358	! cldfra[1-4]D: cloud fraction values [1]
	359	! Ndim: number of dimensions of the input data
	360	! d[1-4]: dimensions of 'cldfra'
	361	! zdim: number of the vertical-dimension within the matrix
	362	! clt1D: total cloudiness for the 1D cldfra
	363	! clt2D1: total cloudiness for the 2D cldfra and d1 being 'zdim'
	364	! clt2D2: total cloudiness for the 2D cldfra and d2 being 'zdim'
	365	! clt3D1: total cloudiness for the 3D cldfra and d1 being 'zdim'
	366	! clt3D2: total cloudiness for the 3D cldfra and d2 being 'zdim'
	367	! clt3D3: total cloudiness for the 3D cldfra and d3 being 'zdim'
	368	! clt4D1: total cloudiness for the 4D cldfra and d1 being 'zdim'
	369	! clt4D2: total cloudiness for the 4D cldfra and d2 being 'zdim'
	370	! clt4D3: total cloudiness for the 4D cldfra and d3 being 'zdim'
	371	! clt4D4: total cloudiness for the 4D cldfra and d4 being 'zdim'
	372
	373	fname = 'compute_clt'
	374
	375	SELECT CASE (Ndim)
	376	CASE (1)
	377	clt1D = var_clt(cldfra1D, d1)
	378	CASE (2)
	379	IF (zdim == 1) THEN
	380	DO i=1, d2
	381	clt2D1(i) = var_clt(cldfra2D(:,i), d1)
	382	END DO
	383	ELSE IF (zdim == 2) THEN
	384	DO i=1, d1
	385	clt2D2(i) = var_clt(cldfra2D(:,i), d2)
	386	END DO
	387	ELSE
	388	PRINT *,TRIM(ErrWarnMsg('err'))
	389	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	390	PRINT *,' accepted values: 1,2'
	391	STOP
	392	END IF
	393	CASE (3)
	394	IF (zdim == 1) THEN
	395	DO i=1, d2
	396	DO j=1, d3
	397	clt3D1(i,j) = var_clt(cldfra3D(:,i,j), d1)
	398	END DO
	399	END DO
	400	ELSE IF (zdim == 2) THEN
	401	DO i=1, d1
	402	DO j=1, d3
	403	clt3D2(i,j) = var_clt(cldfra3D(i,:,j), d2)
	404	END DO
	405	END DO
	406	ELSE IF (zdim == 3) THEN
	407	DO i=1, d1
	408	DO j=1, d2
	409	clt3D3(i,j) = var_clt(cldfra3D(i,j,:), d3)
	410	END DO
	411	END DO
	412	ELSE
	413	PRINT *,TRIM(ErrWarnMsg('err'))
	414	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	415	PRINT *,' accepted values: 1,2,3'
	416	STOP
	417	END IF
	418	CASE (4)
	419	IF (zdim == 1) THEN
	420	DO i=1, d2
	421	DO j=1, d3
	422	DO k=1, d4
	423	clt4D1(i,j,k) = var_clt(cldfra4D(:,i,j,k), d1)
	424	END DO
	425	END DO
	426	END DO
	427	ELSE IF (zdim == 2) THEN
	428	DO i=1, d1
	429	DO j=1, d3
	430	DO k=1, d4
	431	clt4D2(i,j,k) = var_clt(cldfra4D(i,:,j,k), d2)
	432	END DO
	433	END DO
	434	END DO
	435	ELSE IF (zdim == 3) THEN
	436	DO i=1, d2
	437	DO j=1, d3
	438	DO k=1, d4
	439	clt4D3(i,j,k) = var_clt(cldfra4D(i,j,:,k), d3)
	440	END DO
	441	END DO
	442	END DO
	443	ELSE IF (zdim == 4) THEN
	444	DO i=1, d1
	445	DO j=1, d2
	446	DO k=1, d3
	447	clt4D4(i,j,k) = var_clt(cldfra4D(i,j,k,:), d4)
	448	END DO
	449	END DO
	450	END DO
	451	ELSE
	452	PRINT *,TRIM(ErrWarnMsg('err'))
	453	PRINT *,' ' // TRIM(fname) // ': wrong zdim:', zdim,' for Ndim=', Ndim, ' !!'
	454	PRINT *,' accepted values: 1,2,3,4'
	455	STOP
	456	END IF
	457	CASE DEFAULT
	458	PRINT *,TRIM(ErrWarnMsg('err'))
	459	PRINT *,' ' // TRIM(fname) // ': Ndim:', Ndim,' not ready !!'
	460	STOP
	461	END SELECT
	462
	463	RETURN
	464
	465	END SUBROUTINE compute_clt
	466
[1762]	467	SUBROUTINE compute_massvertint1D(var, mutot, dz, deta, integral)
	468	! Subroutine to vertically integrate a 1D variable in eta vertical coordinates
	469
	470	IMPLICIT NONE
	471
	472	INTEGER, INTENT(in) :: dz
	473	REAL(r_k), INTENT(in) :: mutot
	474	REAL(r_k), DIMENSION(dz), INTENT(in) :: var, deta
	475	REAL(r_k), INTENT(out) :: integral
	476
	477	! Local
	478	INTEGER :: k
	479
	480	!!!!!!! Variables
	481	! var: vertical variable to integrate (assuming kgkg-1)
	482	! mutot: total dry-air mass in column
	483	! dz: vertical dimension
	484	! deta: eta-levels difference between full eta-layers
	485
	486	fname = 'compute_massvertint1D'
	487
	488	! integral=0.
	489	! DO k=1,dz
	490	! integral = integral + var(k)*deta(k)
	491	! END DO
	492	integral = SUM(var*deta)
	493
	494	integral=integral*mutot/g
	495
	496	RETURN
	497
	498	END SUBROUTINE compute_massvertint1D
	499
	500	SUBROUTINE compute_zint4D(var4D, dlev, zweight, d1, d2, d3, d4, int3D)
	501	! Subroutine to vertically integrate a 4D variable in any vertical coordinates
	502
	503	IMPLICIT NONE
	504
	505	INTEGER, INTENT(in) :: d1,d2,d3,d4
	506	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: var4D, dlev, zweight
	507	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: int3D
	508
	509	! Local
	510	INTEGER :: i,j,l
	511
	512	!!!!!!! Variables
	513	! var4D: vertical variable to integrate
	514	! dlev: height of layers
	515	! zweight: weight for each level to be applied (=1. for no effect)
	516
	517	fname = 'compute_zint4D'
	518
	519	DO i=1,d1
	520	DO j=1,d2
	521	DO l=1,d4
	522	CALL compute_vertint1D(var4D(i,j,:,l),d3, dlev(i,j,:,l), zweight(i,j,:,l), &
	523	int3D(i,j,l))
	524	END DO
	525	END DO
	526	END DO
	527
	528	RETURN
	529
	530	END SUBROUTINE compute_zint4D
	531
	532	SUBROUTINE compute_vertint1D(var, dz, deta, zweight, integral)
	533	! Subroutine to vertically integrate a 1D variable in any vertical coordinates
	534
	535	IMPLICIT NONE
	536
	537	INTEGER, INTENT(in) :: dz
	538	REAL(r_k), DIMENSION(dz), INTENT(in) :: var, deta, zweight
	539	REAL(r_k), INTENT(out) :: integral
	540
	541	! Local
	542	INTEGER :: k
	543
	544	!!!!!!! Variables
	545	! var: vertical variable to integrate
	546	! dz: vertical dimension
	547	! deta: eta-levels difference between layers
	548	! zweight: weight for each level to be applied (=1. for no effect)
	549
	550	fname = 'compute_vertint1D'
	551
	552	! integral=0.
	553	! DO k=1,dz
	554	! integral = integral + var(k)*deta(k)
	555	! END DO
	556	integral = SUM(vardetazweight)
	557
	558	RETURN
	559
	560	END SUBROUTINE compute_vertint1D
	561
[1759]	562	SUBROUTINE compute_cape_afwa4D(ta, hur, press, zg, hgt, cape, cin, zlfc, plfc, li, parcelmethod, &
	563	d1, d2, d3, d4)
	564	! Subroutine to use WRF phys/module_diag_afwa.F `buyoancy' subroutine to compute CAPE, CIN, ZLFC, PLFC, LI
	565
	566	IMPLICIT NONE
	567
	568	INTEGER, INTENT(in) :: d1, d2, d3, d4, parcelmethod
	569	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: ta, hur, press, zg
	570	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: hgt
	571	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: cape, cin, zlfc, plfc, li
	572
	573	! Local
	574	INTEGER :: i, j, it
	575	INTEGER :: ofunc
	576
	577	!!!!!!! Variables
	578	! ta: air temperature [K]
	579	! hur: relative humidity [%]
	580	! press: air pressure [Pa]
	581	! zg: geopotential height [gpm]
	582	! hgt: topographical height [m]
	583	! cape: Convective available potential energy [Jkg-1]
	584	! cin: Convective inhibition [Jkg-1]
	585	! zlfc: height at the Level of free convection [m]
	586	! plfc: pressure at the Level of free convection [Pa]
	587	! li: lifted index [1]
	588	! parcelmethod:
	589	! Most Unstable = 1 (default)
	590	! Mean layer = 2
	591	! Surface based = 3
	592
	593	fname = 'compute_cape_afwa4D'
	594
	595	DO i=1, d1
	596	DO j=1, d2
	597	DO it=1, d4
	598	ofunc = var_cape_afwa1D(d3, ta(i,j,:,it), hur(i,j,:,it), press(i,j,:,it), zg(i,j,:,it), &
	599	1, cape(i,j,it), cin(i,j,it), zlfc(i,j,it), plfc(i,j,it), li(i,j,it), parcelmethod)
[1833]	600	IF (zlfc(i,j,it) /= -1.) zlfc(i,j,it) = zlfc(i,j,it) - hgt(i,j)
[1759]	601	END DO
	602	END DO
	603	END DO
	604
	605	RETURN
	606
	607	END SUBROUTINE compute_cape_afwa4D
	608
[1795]	609	SUBROUTINE compute_psl_ecmwf(ps, hgt, T, press, unpress, psl, d1, d2, d4)
	610	! Subroutine to compute sea level pressure using ECMWF method following Mats Hamrud and Philippe Courtier [Pa]
	611
	612	IMPLICIT NONE
	613
	614	INTEGER, INTENT(in) :: d1, d2, d4
	615	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(in) :: ps, T, press, unpress
	616	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: hgt
	617	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: psl
	618
	619	! Local
	620	INTEGER :: i, j, it
	621
	622	!!!!!!! Variables
	623	! ps: surface pressure [Pa]
	624	! hgt: terrain height [m]
	625	! T: temperature at first half-mass level [K]
	626	! press: pressure at first full levels [Pa]
	627	! unpress: pressure at first mass (half) levels [Pa]
	628	! psl: sea-level pressure [Pa]
	629
	630	fname = 'compute_psl_ecmwf'
	631
	632	DO i=1, d1
	633	DO j=1, d2
	634	DO it=1, d4
	635	CALL var_psl_ecmwf(ps(i,j,it), hgt(i,j), T(i,j,it), unpress(i,j,it), press(i,j,it), &
	636	psl(i,j,it))
	637	END DO
	638	END DO
	639	END DO
	640
	641	RETURN
	642
	643	END SUBROUTINE compute_psl_ecmwf
	644
[1773]	645	SUBROUTINE compute_zmla_generic4D(tpot, qratio, z, hgt, zmla3D, d1, d2, d3, d4)
	646	! Subroutine to compute pbl-height following a generic method
	647	! from Nielsen-Gammon et al., 2008 J. Appl. Meteor. Clim.
	648	! applied also in Garcia-Diez et al., 2013, QJRMS
	649	! where
	650	! "The technique identifies the ML height as a threshold increase of potential temperature from
	651	! its minimum value within the boundary layer."
	652	! here applied similarly to Garcia-Diez et al. where
	653	! zmla = "...first level where potential temperature exceeds the minimum potential temperature
	654	! reached in the mixed layer by more than 1.5 K"
	655
	656	IMPLICIT NONE
	657
	658	INTEGER, INTENT(in) :: d1, d2, d3, d4
	659	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: tpot, qratio, z
	660	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: hgt
	661	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: zmla3D
	662
	663	! Local
	664	INTEGER :: i, j, it
	665
	666	!!!!!!! Variables
	667	! tpot: potential air temperature [K]
	668	! qratio: water vapour mixing ratio [kgkg-1]
	669	! z: height above sea level [m]
	670	! hgt: terrain height [m]
	671	! zmla3D: boundary layer height from surface [m]
	672
	673	fname = 'compute_zmla_generic4D'
	674
	675	DO i=1, d1
	676	DO j=1, d2
	677	DO it=1, d4
	678	CALL var_zmla_generic(d3, qratio(i,j,:,it), tpot(i,j,:,it), z(i,j,:,it), hgt(i,j), &
	679	zmla3D(i,j,it))
	680	END DO
	681	END DO
	682	END DO
	683
	684	RETURN
	685
	686	END SUBROUTINE compute_zmla_generic4D
	687
[1777]	688	SUBROUTINE compute_zwind4D(ua, va, z, uas, vas, sina, cosa, zextrap, uaz, vaz, d1, d2, d3, d4)
[1776]	689	! Subroutine to compute extrapolate the wind at a given height following the 'power law' methodology
[1773]	690
[1776]	691	IMPLICIT NONE
	692
	693	INTEGER, INTENT(in) :: d1, d2, d3, d4
[1777]	694	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: ua, va, z
[1776]	695	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(in) :: uas, vas
[1777]	696	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: sina, cosa
[1776]	697	REAL(r_k), INTENT(in) :: zextrap
	698	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: uaz, vaz
	699
	700	! Local
	701	INTEGER :: i, j, it
	702
	703	!!!!!!! Variables
	704	! tpot: potential air temperature [K]
	705	! qratio: water vapour mixing ratio [kgkg-1]
[1777]	706	! z: height above surface [m]
[1776]	707	! sina, cosa: local sine and cosine of map rotation [1.]
	708	! zmla3D: boundary layer height from surface [m]
	709
	710	fname = 'compute_zwind4D'
	711
	712	DO i=1, d1
	713	DO j=1, d2
	714	DO it=1, d4
[1777]	715	CALL var_zwind(d3, ua(i,j,:,it), va(i,j,:,it), z(i,j,:,it), uas(i,j,it), vas(i,j,it), &
	716	sina(i,j), cosa(i,j), zextrap, uaz(i,j,it), vaz(i,j,it))
[1776]	717	END DO
	718	END DO
	719	END DO
	720
	721	RETURN
	722
	723	END SUBROUTINE compute_zwind4D
	724
[1784]	725	SUBROUTINE compute_zwind_log4D(ua, va, z, uas, vas, sina, cosa, zextrap, uaz, vaz, d1, d2, d3, d4)
	726	! Subroutine to compute extrapolate the wind at a given height following the 'logarithmic law' methodology
	727
	728	IMPLICIT NONE
	729
	730	INTEGER, INTENT(in) :: d1, d2, d3, d4
	731	REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in) :: ua, va, z
	732	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(in) :: uas, vas
	733	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: sina, cosa
	734	REAL(r_k), INTENT(in) :: zextrap
	735	REAL(r_k), DIMENSION(d1,d2,d4), INTENT(out) :: uaz, vaz
	736
	737	! Local
	738	INTEGER :: i, j, it
	739
	740	!!!!!!! Variables
	741	! tpot: potential air temperature [K]
	742	! qratio: water vapour mixing ratio [kgkg-1]
	743	! z: height above surface [m]
	744	! sina, cosa: local sine and cosine of map rotation [1.]
	745	! zmla3D: boundary layer height from surface [m]
	746
	747	fname = 'compute_zwind_log4D'
	748
	749	DO i=1, d1
	750	DO j=1, d2
	751	DO it=1, d4
	752	CALL var_zwind_log(d3, ua(i,j,:,it), va(i,j,:,it), z(i,j,:,it), uas(i,j,it), vas(i,j,it), &
	753	sina(i,j), cosa(i,j), zextrap, uaz(i,j,it), vaz(i,j,it))
	754	END DO
	755	END DO
	756	END DO
	757
	758	RETURN
	759
	760	END SUBROUTINE compute_zwind_log4D
	761
[1783]	762	SUBROUTINE compute_zwindMO3D(d1, d2, d3, ust, znt, rmol, uas, vas, sina, cosa, newz, uznew, vznew)
	763	! Subroutine to compute extrapolate the wind at a given height following the 'power law' methodology
[1784]	764	! NOTE: only usefull for newz < 80. m
[1783]	765
	766	IMPLICIT NONE
	767
	768	INTEGER, INTENT(in) :: d1, d2, d3
	769	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: ust, znt, rmol
	770	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: uas, vas
	771	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: sina, cosa
	772	REAL(r_k), INTENT(in) :: newz
	773	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out) :: uznew, vznew
	774
	775	! Local
	776	INTEGER :: i, j, it
	777
	778	!!!!!!! Variables
	779	! ust: u* in similarity theory [ms-1]
	780	! znt: thermal time-varying roughness length [m]
	781	! rmol: Inverse of the Obukhov length [m-1]
	782	! uas: x-component 10-m wind speed [ms-1]
	783	! vas: y-component 10-m wind speed [ms-1]
	784	! sina, cosa: local sine and cosine of map rotation [1.]
	785
	786	fname = 'compute_zwindMO3D'
	787
	788	DO i=1, d1
	789	DO j=1, d2
	790	DO it=1, d3
	791	CALL var_zwind_MOtheor(ust(i,j,it), znt(i,j,it), rmol(i,j,it), uas(i,j,it), vas(i,j,it), &
	792	sina(i,j), cosa(i,j), newz, uznew(i,j,it), vznew(i,j,it))
	793	END DO
	794	END DO
	795	END DO
	796
	797	RETURN
	798
	799	END SUBROUTINE compute_zwindMO3D
	800
[1804]	801	SUBROUTINE compute_potevap_orPM3D(d1, d2, d3, rho1, ust, uas, vas, tas, ps, qv1, potevap)
	802	! Subroutine to compute potential evapotranspiration Penman-Monteith formulation implemented in
[1833]	803	! ORCHIDEE in src_sechiba/enerbil.f90
[1804]	804
	805	IMPLICIT NONE
	806
	807	INTEGER, INTENT(in) :: d1, d2, d3
	808	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: rho1, ust, uas, vas, tas, ps, qv1
	809	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out) :: potevap
	810
	811	! Local
	812	INTEGER :: i, j, it
	813
	814	!!!!!!! Variables
	815	! rho1: atsmophere density at the first layer [kgm-3]
	816	! ust: u* in similarity theory [ms-1]
	817	! uas: x-component 10-m wind speed [ms-1]
	818	! vas: y-component 10-m wind speed [ms-1]
	819	! tas: 2-m atmosphere temperature [K]
	820	! ps: surface pressure [Pa]
	821	! qv1: 1st layer atmospheric mixing ratio [kgkg-1]
	822	! potevap: potential evapo transpiration [kgm-2s-1]
	823
	824	fname = 'compute_potevap_orPM3D'
	825
	826	DO i=1, d1
	827	DO j=1, d2
	828	DO it=1, d3
	829	CALL var_potevap_orPM(rho1(i,j,it), ust(i,j,it), uas(i,j,it), vas(i,j,it), tas(i,j,it), &
	830	ps(i,j,it), qv1(i,j,it), potevap(i,j,it))
	831	END DO
	832	END DO
	833	END DO
	834
	835	RETURN
	836
	837	END SUBROUTINE compute_potevap_orPM3D
	838
[1908]	839	SUBROUTINE compute_fog_K84(d1, d2, d3, qc, qi, fog, vis)
	840	! Subroutine to compute fog: qcloud + qice /= 0.
	841	! And visibility following Kunkel, B. A., (1984): Parameterization of droplet terminal velocity and
	842	! extinction coefficient in fog models. J. Climate Appl. Meteor., 23, 34â41.
	843
	844	IMPLICIT NONE
	845
	846	INTEGER, INTENT(in) :: d1, d2, d3
	847	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: qc, qi
	848	INTEGER, DIMENSION(d1,d2,d3), INTENT(out) :: fog
	849	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out) :: vis
	850
	851	! Local
	852	INTEGER :: i, j, it
	853
	854	!!!!!!! Variables
	855	! qc: cloud mixing ratio [kgkg-1]
	856	! qi, ice mixing ratio [kgkg-1]
	857	! fog: presence of fog (1: yes, 0: no)
	858	! vis: visibility within fog [km]
	859
	860	fname = 'compute_fog_K84'
	861
	862	DO i=1, d1
	863	DO j=1, d2
	864	DO it=1, d3
	865	CALL var_fog_K84(qc(i,j,it), qi(i,j,it), fog(i,j,it), vis(i,j,it))
	866	END DO
	867	END DO
	868	END DO
	869
	870	RETURN
	871
	872	END SUBROUTINE compute_fog_K84
	873
[1909]	874	SUBROUTINE compute_fog_RUC(d1, d2, d3, qv, ta, pres, fog, vis)
[1908]	875	! Subroutine to compute fog: qcloud + qice /= 0.
	876	! And visibility following RUC method Smirnova, T. G., S. G. Benjamin, and J. M. Brown, 2000: Case
	877	! study verification of RUC/MAPS fog and visibility forecasts. Preprints, 9 th Conference on
	878	! Aviation, Range, and Aerospace Meteorlogy, AMS, Orlando, FL, Sep. 2000. Paper#2.3, 6 pp.
	879
	880	IMPLICIT NONE
	881
	882	INTEGER, INTENT(in) :: d1, d2, d3
[1909]	883	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: qv, ta, pres
[1908]	884	INTEGER, DIMENSION(d1,d2,d3), INTENT(out) :: fog
	885	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out) :: vis
	886
	887	! Local
	888	INTEGER :: i, j, it
	889
	890	!!!!!!! Variables
[1909]	891	! qv: water vapor mixing ratio [kgkg-1]
	892	! ta: temperature [K]
	893	! pres: pressure [Pa]
[1908]	894	! fog: presence of fog (1: yes, 0: no)
	895	! vis: visibility within fog [km]
	896
[1909]	897	fname = 'compute_fog_RUC'
[1908]	898
	899	DO i=1, d1
	900	DO j=1, d2
	901	DO it=1, d3
[1909]	902	CALL var_fog_RUC(qv(i,j,it), ta(i,j,it), pres(i,j,it), fog(i,j,it), vis(i,j,it))
[1908]	903	END DO
	904	END DO
	905	END DO
	906
	907	RETURN
	908
	909	END SUBROUTINE compute_fog_RUC
	910
[1909]	911	SUBROUTINE compute_fog_FRAML50(d1, d2, d3, qv, ta, pres, fog, vis)
	912	! Subroutine to compute fog (vis < 1 km) and visibility following
	913	! Gultepe, I. and J.A. Milbrandt, 2010: Probabilistic Parameterizations of Visibility Using
	914	! Observations of Rain Precipitation Rate, Relative Humidity, and Visibility. J. Appl. Meteor.
	915	! Climatol., 49, 36-46, https://doi.org/10.1175/2009JAMC1927.1
	916	! Interest is focused on a 'general' fog/visibilty approach, thus the fit at 50 % of probability is
	917	! chosen
	918	! Effects from precipitation are not considered
	919
	920	IMPLICIT NONE
	921
	922	INTEGER, INTENT(in) :: d1, d2, d3
	923	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in) :: qv, ta, pres
	924	INTEGER, DIMENSION(d1,d2,d3), INTENT(out) :: fog
	925	REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out) :: vis
	926
	927	! Local
	928	INTEGER :: i, j, it
	929
	930	!!!!!!! Variables
	931	! qv: mixing ratio in [kgkg-1]
	932	! ta: temperature [K]
	933	! pres: pressure field [Pa]
	934	! fog: presence of fog (1: yes, 0: no)
	935	! vis: visibility within fog [km]
	936
	937	fname = 'compute_fog_FRAML50'
	938
	939	DO i=1, d1
	940	DO j=1, d2
	941	DO it=1, d3
	942	CALL var_fog_FRAML50(qv(i,j,it), ta(i,j,it), pres(i,j,it), fog(i,j,it), vis(i,j,it))
	943	END DO
	944	END DO
	945	END DO
	946
	947	RETURN
	948
	949	END SUBROUTINE compute_fog_FRAML50
	950
[2208]	951	SUBROUTINE compute_range_faces(d1, d2, lon, lat, hgt, face, faces)
	952	! Subroutine to compute faces [uphill, valleys, downhill] of a mountain range along a given face
	953
	954	IMPLICIT NONE
	955
	956	INTEGER, INTENT(in) :: d1, d2
	957	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: lon, lat, hgt
	958	CHARACTER(len=*) :: face
	959	INTEGER, DIMENSION(d1,d2), INTENT(out) :: faces
	960
	961	! Local
	962	INTEGER :: i, j
	963
	964	!!!!!!! Variables
	965	! lon: longitude [degrees east]
	966	! lat: latitude [degrees north]
	967	! hgt: topograpical height [m]
	968
	969	fname = 'compute_range_faces'
	970
	971	IF (TRIM(face) == 'WE') THEN
	972	DO j=1, d2
	973	CALL var_range_faces(d1, lon(:,j), lat(:,j), hgt(:,j), faces(:,j))
	974	END DO
	975	ELSE IF (TRIM(face) == 'SN') THEN
	976	DO i=1, d1
	977	CALL var_range_faces(d2, lon(i,:), lat(i,:), hgt(i,:), faces(i,:))
	978	END DO
	979	ELSE
	980	PRINT *,TRIM(ErrWarnMsg('err'))
	981	PRINT *,' ' // TRIM(fname) // ": wrong face: '" // TRIM(face) // "' !!"
	982	PRINT *,' accepted ones: WE, SN'
	983	STOP
	984	END IF
	985
	986	RETURN
	987
	988	END SUBROUTINE compute_range_faces
	989
[770]	990	END MODULE module_ForDiagnostics

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format