SURFEX v8.1
General documentation of Surfex
facdec.F90
Go to the documentation of this file.
1 ! Oct-2012 P. Marguinaud 64b LFI
2 ! Jan-2011 P. Marguinaud Thread-safe FA
3 SUBROUTINE facdec_fort &
4 & (fa, krep, pa, pmin, knbit, kdec)
5 USE fa_mod, ONLY : fa_com
6 USE parkind1, ONLY : jprb
7 USE yomhook , ONLY : lhook, dr_hook
9 IMPLICIT NONE
10 !****
11 ! Sous-programme de calcul du FACteur d'echelle DECimal optimal
12 ! pour coder un champ d'amplitude PA donnee, sur KNBIT bits avec
13 ! une precision maximale, compte tenu d'un intervalle fixe
14 ! de facteurs decimaux.
15 !**
16 ! Arguments : KREP (Sortie) ==> Code-reponse du sous-programme;
17 ! 0: OK; -1: pb rencontre
18 ! PA (Entree) ==> Amplitude du champ a compacter;
19 ! PMIN (Entree) ==> Minimum du champ a compacter;
20 ! KNBIT (Entree) ==> Nb de bits servant au compactage;
21 ! KDEC (Sortie) ==> Facteur d'echelle decimal optimal;
22 !
23 ! Modifications:
24 ! R. El Khatib 12-Aug-2009 Bugfix for compatibility with Gribex
25 ! R. El Khatib 22-Jul-2015 Threshold for the field amplitude raised
26 ! from 1.E-12 to 1.E-11 in order to avoid the deterioration of too flat fields
27 !
28 TYPE(fa_com) :: FA
29 REAL (KIND=JPDBLR) PA, PMIN
30 INTEGER (KIND=JPLIKB) KREP, KNBIT, KDEC
31 !
32 REAL (KIND=JPDBLR) XNBINT, XTINYR4, XHUGER4
33 INTEGER (KIND=JPLIKB) IDECMIN, IDECMAX, INBINT
34 INTEGER (KIND=JPLIKB) JDEC, IE, INUTIL, INUMAX, IEMAX
35 INTEGER (KIND=JPLIKB) INU0, IE0
36 
37 !
38 !**
39 ! 1. - CONTROLES ET INITIALISATIONS.
40 !-----------------------------------------------------------------------
41 !
42 REAL(KIND=JPRB) :: ZHOOK_HANDLE
43 IF (lhook) CALL dr_hook('FACDEC_MT',0,zhook_handle)
44 IF (knbit.LE.0 .OR. knbit.GT.64) THEN
45  krep=-1
46  WRITE (unit=fa%NULOUT,fmt=*)'****'
47  WRITE (unit=fa%NULOUT,fmt=*) &
48 & '**** FACDEC: ERROR, bits number out of range 1-64'
49  WRITE (unit=fa%NULOUT,fmt=*)'**** KNBIT = ',knbit
50  WRITE (unit=fa%NULOUT,fmt=*) &
51 & '**** ! Optimal decimal scale factor is not computed !'
52  WRITE (unit=fa%NULOUT,fmt=*)'****'
53  GOTO 1001
54 ENDIF
55 IF (abs(pa) .LT. tiny(pa)) THEN
56  krep = 0
57  kdec = 0
58  IF (fa%LFAMOP) THEN
59  WRITE (unit=fa%NULOUT,fmt=*)'////'
60  WRITE (unit=fa%NULOUT,fmt=*) &
61 & '//// FACDEC: WARNING, range of the field is null :', &
62 & pa
63  WRITE (unit=fa%NULOUT,fmt=*)'////'
64  ENDIF
65  GOTO 1001
66 ENDIF
67 
68 krep = 0
69 inumax = 0
70 idecmin = -15
71 idecmax = 5
72 xtinyr4 = tiny(1._4)
73 xhuger4 = huge(1._4)
74 inbint = 2**knbit -1
75 xnbint = real(inbint, jpdblr)
76 ! Cas du facteur decimal nul (reference a calculer dans tous les cas)
77 jdec = 0
78 CALL factec_fort &
79 & (fa, krep, pa, knbit, jdec, ie0, inu0)
80 !
81 !**
82 ! 2. - BOUCLE TESTANT LES FACTEURS DECIMAUX
83 !-----------------------------------------------------------------------
84 !
85 !!OCL SCALAR
86 DO jdec = idecmin, idecmax
87 !
88 ! 3 tests sur la pertinence de ce facteur decimal
89 !
90 ! 0/ PA * 10**JDEC > 1.E-11
91 ! Du fait d'un test exagere dans gribex
92  IF (pa * 10._jpdblr**jdec .LE. 1.e-11_jpdblr) cycle
93 ! 1/ PMIN * 10**JDEC > TINY(real*4) pour decodage eventuel
94 ! de la valeur de reference en arithmetique 32 bits IEEE.
95  IF (abs(pmin) .GT. tiny(pmin)) THEN
96  IF (log10(abs(pmin))+REAL (JDEC, JPDBLR) .LE. LOG10( xtinyr4 ) ) cycle
97  ELSE
98 ! instruction "bidon" pour que l'optimisation ne fasse pas
99 ! l'evaluation de LOG10(ABS(PMIN)) par anticipation (Pb sur VPP5000)
100  pmin = 0._jpdblr
101  ENDIF
102 ! 2/ PA * 10**JDEC inclus dans le domaine de validite des real*8
103  IF (abs(log10(abs(pa))+real(jdec, jpdblr)) .GE. REAL (RANGE(PA), JPDBLR)) cycle
104 !
105  CALL factec_fort &
106 & (fa, krep, pa, knbit, jdec, ie, inutil)
107  IF (krep.NE.0) cycle
108 ! 3/ PMIN*10**JDEC + (2**KNBIT-1)*2**IE < HUGE(real*4) pour decodage
109 ! eventuel du max du champ en arithmetique 32 bits IEEE.
110  IF (pmin*10._jpdblr**jdec + xnbint*2._jpdblr**ie .GE. xhuger4) cycle
111 ! 4/ Le facteur d'echelle binaire doit pouvoir etre code sur 1 octet,
112 ! il doit donc etre compris entre -126 et 127:
113  IF (ie.LT.-126 .OR. ie.GT.127) cycle
114 !
115  IF (inutil.GT.inumax) THEN
116  inumax = inutil
117  iemax = ie
118  kdec = jdec
119  ENDIF
120 ENDDO
121 !
122 IF (inumax.EQ.0) THEN
123  krep=-1
124  WRITE (unit=fa%NULOUT,fmt=*)'****'
125  WRITE (unit=fa%NULOUT,fmt=*) &
126 & '**** FACDEC: all the decimal factors comprised between'
127  WRITE (unit=fa%NULOUT,fmt=*) &
128 & '**** ',idecmin,' and ',idecmax,' are rejected !!'
129  WRITE (unit=fa%NULOUT,fmt=*) &
130 & '**** Range and min of the field are :', pa, pmin
131  WRITE (unit=fa%NULOUT,fmt=*) &
132 & '****'
133  WRITE (unit=fa%NULOUT,fmt=*)'****'
134  GOTO 1001
135 ENDIF
136 IF (fa%LFAMOP) THEN
137  WRITE (unit=fa%NULOUT,fmt=*) &
138 & 'FACDEC: champ d''amplitude ',pa,' ,de minimum ',pmin
139  WRITE (unit=fa%NULOUT,fmt=*) &
140 & ' => fact decimal opt de',kdec, &
141 & ' ,pour 1 fact binaire de ',iemax
142  WRITE (unit=fa%NULOUT,fmt='(1X,A,I3,A,I9,A,I9,A,F5.1,A,E11.4)') &
143 & ' Eff des', &
144 & knbit,' bits = ',inumax,' sur ',inbint,' soit: ', &
145 & real(inumax, jpdblr)*100._jpdblr/xnbint, &
146 & ' % et une precision de ', &
147 & 10._jpdblr**(-kdec)*2._jpdblr**(iemax-1)
148  WRITE (unit=fa%NULOUT,fmt=*) &
149 & ' a comparer, si le fact decimal = 0, avec'
150  WRITE (unit=fa%NULOUT,fmt='(1X,A,I9,A,I9,A,F5.1,A,E11.4)') &
151 & ' une efficacite de ', &
152 & inu0,' sur ',inbint,' soit: ', &
153 & real(inu0, jpdblr)*100._jpdblr/xnbint, &
154 & ' % et une precision de ',2._jpdblr**(ie0-1)
155 ENDIF
156 !
157 1001 CONTINUE
158 !
159 IF (lhook) CALL dr_hook('FACDEC_MT',1,zhook_handle)
160 END SUBROUTINE facdec_fort
161 
162 
163 
164 ! Oct-2012 P. Marguinaud 64b LFI
165 SUBROUTINE facdec64 &
166 & (krep, pa, pmin, knbit, kdec)
167 USE fa_mod, ONLY : fa => fa_com_default, &
170 USE lfi_precision
171 IMPLICIT NONE
172 ! Arguments
173 INTEGER (KIND=JPLIKB) KREP ! OUT
174 REAL (KIND=JPDBLR) PA ! IN
175 REAL (KIND=JPDBLR) PMIN ! IN
176 INTEGER (KIND=JPLIKB) KNBIT ! IN
177 INTEGER (KIND=JPLIKB) KDEC ! OUT
178 
179 IF (.NOT. fa_com_default_init) CALL new_fa_default ()
180 
181 CALL facdec_fort &
182 & (fa, krep, pa, pmin, knbit, kdec)
183 
184 END SUBROUTINE facdec64
185 
186 SUBROUTINE facdec &
187 & (krep, pa, pmin, knbit, kdec)
188 USE fa_mod, ONLY : fa => fa_com_default, &
191 USE lfi_precision
192 IMPLICIT NONE
193 ! Arguments
194 INTEGER (KIND=JPLIKM) KREP ! OUT
195 REAL (KIND=JPDBLR) PA ! IN
196 REAL (KIND=JPDBLR) PMIN ! IN
197 INTEGER (KIND=JPLIKM) KNBIT ! IN
198 INTEGER (KIND=JPLIKM) KDEC ! OUT
199 
200 IF (.NOT. fa_com_default_init) CALL new_fa_default ()
201 
202 CALL facdec_mt &
203 & (fa, krep, pa, pmin, knbit, kdec)
204 
205 END SUBROUTINE facdec
206 
207 SUBROUTINE facdec_mt &
208 & (fa, krep, pa, pmin, knbit, kdec)
209 USE fa_mod, ONLY : fa_com
210 USE lfi_precision
211 IMPLICIT NONE
212 ! Arguments
213 type(fa_com) fa ! INOUT
214 INTEGER (KIND=JPLIKM) KREP ! OUT
215 REAL (KIND=JPDBLR) PA ! IN
216 REAL (KIND=JPDBLR) PMIN ! IN
217 INTEGER (KIND=JPLIKM) KNBIT ! IN
218 INTEGER (KIND=JPLIKM) KDEC ! OUT
219 ! Local integers
220 INTEGER (KIND=JPLIKB) IREP ! OUT
221 INTEGER (KIND=JPLIKB) INBIT ! IN
222 INTEGER (KIND=JPLIKB) IDEC ! OUT
223 ! Convert arguments
224 
225 inbit = int( knbit, jplikb)
226 
227 CALL facdec_fort &
228 & (fa, irep, pa, pmin, inbit, idec)
229 
230 krep = int( irep, jplikm)
231 kdec = int( idec, jplikm)
232 
233 END SUBROUTINE facdec_mt
234 
235 !INTF KREP OUT
236 !INTF PA IN
237 !INTF PMIN IN
238 !INTF KNBIT IN
239 !INTF KDEC OUT
subroutine factec_fort(FA, KREP, PA, KNBIT, KDEC, KE, KNUTIL)
Definition: factec.F90:5
subroutine facdec(KREP, PA, PMIN, KNBIT, KDEC)
Definition: facdec.F90:188
integer, parameter jplikb
subroutine facdec_fort(FA, KREP, PA, PMIN, KNBIT, KDEC)
Definition: facdec.F90:5
subroutine facdec64(KREP, PA, PMIN, KNBIT, KDEC)
Definition: facdec.F90:167
logical, save fa_com_default_init
Definition: fa_mod.F90:477
subroutine new_fa_default()
Definition: fa_mod.F90:649
subroutine facdec_mt(FA, KREP, PA, PMIN, KNBIT, KDEC)
Definition: facdec.F90:209
Definition: fa_mod.F90:1
integer, parameter jprb
Definition: parkind1.F90:32
integer, parameter jpdblr
logical lhook
Definition: yomhook.F90:15
integer, parameter jplikm
type(fa_com), target, save fa_com_default
Definition: fa_mod.F90:476
real8 real
Definition: privpub.h:396