XBeach
|
00001 module means_module 00002 00003 ! This module uses the On-line algorithm according to Knuth (1998) to 00004 ! determine the variance on the fly 00005 ! Note: this approximation can be improved. 00006 use xmpi_module 00007 use mnemmodule, only: arraytype, maxnamelen, chartoindex 00008 00009 implicit none 00010 save 00011 00012 private 00013 public meanspars, meansparsglobal, meansparslocal, means_init, makeaverage, clearaverage 00014 public makecrossvector 00015 #ifdef USEMPI 00016 public means_collect 00017 #endif 00018 type meanspars 00019 ! Name 00020 character(maxnamelen) :: name 00021 ! Rank 00022 integer :: rank 00023 00024 ! Array type 00025 type(arraytype) :: t 00026 00027 ! Keep time average variables 00028 real*8,dimension(:,:), pointer:: mean2d => NULL() 00029 real*8,dimension(:,:,:), pointer:: mean3d => NULL() 00030 real*8,dimension(:,:,:,:), pointer:: mean4d=> NULL() 00031 ! Keep variance of variables 00032 real*8,dimension(:,:), pointer:: variance2d => NULL() 00033 real*8,dimension(:,:,:), pointer:: variance3d => NULL() 00034 real*8,dimension(:,:,:,:), pointer:: variance4d => NULL() 00035 ! Needed for variance calculation 00036 real*8,dimension(:,:), pointer:: variancecrossterm2d => NULL() 00037 real*8,dimension(:,:,:), pointer:: variancecrossterm3d => NULL() 00038 real*8,dimension(:,:,:,:), pointer:: variancecrossterm4d=> NULL() 00039 ! Needed for variance calculation 00040 real*8,dimension(:,:), pointer:: variancesquareterm2d=> NULL() 00041 real*8,dimension(:,:,:), pointer:: variancesquareterm3d => NULL() 00042 real*8,dimension(:,:,:,:), pointer:: variancesquareterm4d => NULL() 00043 ! Keep time min variables 00044 real*8,dimension(:,:), pointer:: min2d => NULL() 00045 real*8,dimension(:,:,:), pointer:: min3d => NULL() 00046 real*8,dimension(:,:,:,:), pointer:: min4d => NULL() 00047 ! Keep time max variables 00048 real*8,dimension(:,:), pointer:: max2d => NULL() 00049 real*8,dimension(:,:,:), pointer:: max3d => NULL() 00050 real*8,dimension(:,:,:,:), pointer:: max4d => NULL() 00051 end type meanspars 00052 00053 type(meanspars),dimension(:),allocatable :: meansparsglobal 00054 type(meanspars),dimension(:),allocatable :: meansparslocal 00055 00056 contains 00057 00058 subroutine means_init(sg,sl,par) 00059 use spaceparams 00060 use params, only: parameters 00061 use mnemmodule 00062 00063 implicit none 00064 00065 type(spacepars),intent(inout) :: sg 00066 type(spacepars),intent(in) :: sl 00067 type(parameters),intent(in) :: par 00068 type(arraytype) :: t 00069 integer :: i,index,d1,d2,d3,d4 00070 00071 d1 = -123 00072 d2 = -123 00073 d3 = -123 00074 d4 = -123 00075 00076 if (par%nmeanvar>0) then 00077 00078 allocate(meansparsglobal(par%nmeanvar)) 00079 allocate(meansparslocal(par%nmeanvar)) 00080 00081 do i=1,par%nmeanvar 00082 index=chartoindex(par%meanvars(i)) 00083 #ifdef USEMPI 00084 if(xomaster) then 00085 ! just to make sure that d1 .. d4 get the correct values: 00086 call index_allocate(sg,par,index,'r') 00087 endif 00088 #endif 00089 call indextos(sg,index,t) 00090 meansparsglobal(i)%name=t%name 00091 meansparsglobal(i)%rank=t%rank 00092 meansparsglobal(i)%t = t 00093 select case (t%rank) 00094 case (2) 00095 if (t%type == 'r') then 00096 #ifdef USEMPI 00097 !t%r2 = 0 ! wwvv todo why? 00098 #endif 00099 d1 = size(t%r2,1) 00100 d2 = size(t%r2,2) 00101 elseif (t%type == 'i') then 00102 #ifdef USEMPI 00103 !t%i2 = 0 00104 #endif 00105 d1 = size(t%i2,1) 00106 d2 = size(t%i2,2) 00107 else 00108 call halt_program 00109 end if 00110 allocate(meansparsglobal(i)%mean2d(d1,d2)) 00111 allocate(meansparsglobal(i)%variance2d(d1,d2)) 00112 allocate(meansparsglobal(i)%variancecrossterm2d(d1,d2)) 00113 allocate(meansparsglobal(i)%variancesquareterm2d(d1,d2)) 00114 allocate(meansparsglobal(i)%min2d(d1,d2)) 00115 allocate(meansparsglobal(i)%max2d(d1,d2)) 00116 meansparsglobal(i)%mean2d = 0.d0 00117 meansparsglobal(i)%variance2d = 0.d0 00118 meansparsglobal(i)%variancecrossterm2d = 0.d0 00119 meansparsglobal(i)%variancesquareterm2d = 0.d0 00120 meansparsglobal(i)%min2d = huge(0.d0) 00121 meansparsglobal(i)%max2d = -1.d0*huge(0.d0) 00122 case (3) 00123 if (t%type == 'r') then 00124 #ifdef USEMPI 00125 t%r3 = 0 00126 #endif 00127 d1 = size(t%r3,1) 00128 d2 = size(t%r3,2) 00129 d3 = size(t%r3,3) 00130 elseif (t%type == 'i') then 00131 #ifdef USEMPI 00132 t%i3 = 0 00133 #endif 00134 d1 = size(t%i3,1) 00135 d2 = size(t%i3,2) 00136 d3 = size(t%i3,3) 00137 else 00138 call halt_program 00139 end if 00140 allocate(meansparsglobal(i)%mean3d(d1,d2,d3)) 00141 allocate(meansparsglobal(i)%variance3d(d1,d2,d3)) 00142 allocate(meansparsglobal(i)%variancecrossterm3d(d1,d2,d3)) 00143 allocate(meansparsglobal(i)%variancesquareterm3d(d1,d2,d3)) 00144 allocate(meansparsglobal(i)%min3d(d1,d2,d3)) 00145 allocate(meansparsglobal(i)%max3d(d1,d2,d3)) 00146 meansparsglobal(i)%mean3d = 0.d0 00147 meansparsglobal(i)%variance3d = 0.d0 00148 meansparsglobal(i)%variancecrossterm3d = 0.d0 00149 meansparsglobal(i)%variancesquareterm3d = 0.d0 00150 meansparsglobal(i)%min3d = huge(0.d0) 00151 meansparsglobal(i)%max3d = -1.d0*huge(0.d0) 00152 case (4) 00153 if (t%type == 'r') then 00154 #ifdef USEMPI 00155 t%r4 = 0 00156 #endif 00157 d1 = size(t%r4,1) 00158 d2 = size(t%r4,2) 00159 d3 = size(t%r4,3) 00160 d4 = size(t%r4,4) 00161 elseif (t%type == 'i') then 00162 #ifdef USEMPI 00163 t%i4 = 0 00164 #endif 00165 d1 = size(t%i4,1) 00166 d2 = size(t%i4,2) 00167 d3 = size(t%i4,3) 00168 d4 = size(t%i4,4) 00169 else 00170 call halt_program 00171 end if 00172 allocate(meansparsglobal(i)%mean4d(d1,d2,d3,d4)) 00173 allocate(meansparsglobal(i)%variance4d(d1,d2,d3,d4)) 00174 allocate(meansparsglobal(i)%variancecrossterm4d(d1,d2,d3,d4)) 00175 allocate(meansparsglobal(i)%variancesquareterm4d(d1,d2,d3,d4)) 00176 allocate(meansparsglobal(i)%min4d(d1,d2,d3,d4)) 00177 allocate(meansparsglobal(i)%max4d(d1,d2,d3,d4)) 00178 meansparsglobal(i)%mean4d = 0.d0 00179 meansparsglobal(i)%variance4d = 0.d0 00180 meansparsglobal(i)%variancecrossterm4d = 0.d0 00181 meansparsglobal(i)%variancesquareterm4d = 0.d0 00182 meansparsglobal(i)%min4d = huge(0.d0) 00183 meansparsglobal(i)%max4d = -1.d0*huge(0.d0) 00184 end select 00185 #ifdef USEMPI 00186 call indextos(sl,index,t) 00187 #else 00188 ! do nothing, we use s global here, no local s available 00189 #endif 00190 meansparslocal(i)%name=t%name 00191 meansparslocal(i)%rank=t%rank 00192 meansparslocal(i)%t = t 00193 select case (t%rank) 00194 case (2) 00195 if (t%type == 'r') then 00196 d1 = size(t%r2,1) 00197 d2 = size(t%r2,2) 00198 elseif (t%type == 'i') then 00199 d1 = size(t%r2,1) 00200 d2 = size(t%r2,2) 00201 else 00202 call halt_program 00203 end if 00204 allocate(meansparslocal(i)%mean2d(d1,d2)) 00205 allocate(meansparslocal(i)%variance2d(d1,d2)) 00206 allocate(meansparslocal(i)%variancecrossterm2d(d1,d2)) 00207 allocate(meansparslocal(i)%variancesquareterm2d(d1,d2)) 00208 allocate(meansparslocal(i)%min2d(d1,d2)) 00209 allocate(meansparslocal(i)%max2d(d1,d2)) 00210 meansparslocal(i)%mean2d = 0.d0 00211 meansparslocal(i)%variance2d = 0.d0 00212 meansparslocal(i)%variancecrossterm2d = 0.d0 00213 meansparslocal(i)%variancesquareterm2d = 0.d0 00214 meansparslocal(i)%min2d = huge(0.d0) 00215 meansparslocal(i)%max2d = -1.d0*huge(0.d0) 00216 case (3) 00217 if (t%type == 'r') then 00218 d1 = size(t%r3,1) 00219 d2 = size(t%r3,2) 00220 d3 = size(t%r3,3) 00221 elseif (t%type == 'i') then 00222 d1 = size(t%r3,1) 00223 d2 = size(t%r3,2) 00224 d3 = size(t%r3,3) 00225 else 00226 call halt_program 00227 end if 00228 allocate(meansparslocal(i)%mean3d(d1,d2,d3)) 00229 allocate(meansparslocal(i)%variance3d(d1,d2,d3)) 00230 allocate(meansparslocal(i)%variancecrossterm3d(d1,d2,d3)) 00231 allocate(meansparslocal(i)%variancesquareterm3d(d1,d2,d3)) 00232 allocate(meansparslocal(i)%min3d(d1,d2,d3)) 00233 allocate(meansparslocal(i)%max3d(d1,d2,d3)) 00234 meansparslocal(i)%mean3d = 0.d0 00235 meansparslocal(i)%variance3d = 0.d0 00236 meansparslocal(i)%variancecrossterm3d = 0.d0 00237 meansparslocal(i)%variancesquareterm3d = 0.d0 00238 meansparslocal(i)%min3d = huge(0.d0) 00239 meansparslocal(i)%max3d = -1.d0*huge(0.d0) 00240 case (4) 00241 if (t%type == 'r') then 00242 d1 = size(t%r4,1) 00243 d2 = size(t%r4,2) 00244 d3 = size(t%r4,3) 00245 d4 = size(t%r4,4) 00246 elseif (t%type == 'i') then 00247 d1 = size(t%r4,1) 00248 d2 = size(t%r4,2) 00249 d3 = size(t%r4,3) 00250 d4 = size(t%r4,4) 00251 else 00252 call halt_program 00253 end if 00254 allocate(meansparslocal(i)%mean4d(d1,d2,d3,d4)) 00255 allocate(meansparslocal(i)%variance4d(d1,d2,d3,d4)) 00256 allocate(meansparslocal(i)%variancecrossterm4d(d1,d2,d3,d4)) 00257 allocate(meansparslocal(i)%variancesquareterm4d(d1,d2,d3,d4)) 00258 allocate(meansparslocal(i)%min4d(d1,d2,d3,d4)) 00259 allocate(meansparslocal(i)%max4d(d1,d2,d3,d4)) 00260 meansparslocal(i)%mean4d = 0.d0 00261 meansparslocal(i)%variance4d = 0.d0 00262 meansparslocal(i)%variancecrossterm4d = 0.d0 00263 meansparslocal(i)%variancesquareterm4d = 0.d0 00264 meansparslocal(i)%min4d = huge(0.d0) 00265 meansparslocal(i)%max4d = -1.d0*huge(0.d0) 00266 end select 00267 00268 enddo 00269 endif 00270 ! wwvv to avoid warning about unused sl: 00271 if (sl%nx .eq. -1) return 00272 end subroutine means_init 00273 00274 00275 subroutine makeaverage(sl,par) 00276 00277 use params, only: parameters 00278 use spaceparams 00279 use mnemmodule, only: chartoindex 00280 use logging_module 00281 use postprocessmod, only: gridrotate 00282 use timestep_module 00283 00284 IMPLICIT NONE 00285 00286 type(parameters), intent(in) :: par 00287 type(spacepars), intent(inout) :: sl 00288 00289 ! keep track of which mean variables are used 00290 integer :: index 00291 integer :: i 00292 real*8 :: mult 00293 00294 type(arraytype) :: t 00295 00296 integer,dimension(sl%nx+1,sl%ny+1) :: tvar2di 00297 real*8,dimension(:,:),allocatable,save :: tvar2d_sin,tvar2d_cos 00298 ! real*8,dimension(sl%nx+1,sl%ny+1) :: 00299 00300 real*8, dimension(sl%nx+1,sl%ny+1) :: oldmean2d,tvar2d 00301 real*8, dimension(:,:,:), allocatable :: oldmean3d,tvar3d 00302 real*8, dimension(:,:,:,:), allocatable :: oldmean4d,tvar4d 00303 00304 real*8, parameter :: numeps = epsilon(0.d0) 00305 logical,save :: initialisedtvarsin = .false. 00306 00307 00308 00309 ! avgtime = tpar%tpm(tpar%itm+1)-tpar%tpm(tpar%itm) 00310 ! updated in varoutput, not needed to calculate here 00311 mult = max(par%dt/par%tintm,0.d0) ! Catch initialization at t=0 00312 00313 do i=1,par%nmeanvar 00314 index=chartoindex(par%meanvars(i)) 00315 call indextos(sl,index,t) 00316 select case (t%rank) 00317 case (2) 00318 oldmean2d=meansparslocal(i)%mean2d 00319 ! Robert: One where, elsewhere, endwhere statement leads to memory leak. 00320 ! maybe ifort bug ? 00321 ! wwvv 2014-10-21 changed tiny into epsilon in this subroutine 00322 where (oldmean2d<numeps .and. oldmean2d>=0.d0) 00323 oldmean2d=numeps 00324 endwhere 00325 where (oldmean2d>-numeps .and. oldmean2d<0.d0) 00326 oldmean2d=-numeps 00327 endwhere 00328 ! Some variables (vectors) are rotated to N-S and E-W direction 00329 if (t%type=='i') then 00330 call gridrotate(t,tvar2di) ! wwvv-todo 00331 tvar2d=dble(tvar2di) 00332 else 00333 call gridrotate(par, sl,t,tvar2d) 00334 endif 00335 if (par%meanvars(i)=='thetamean') then 00336 if (.not. initialisedtvarsin) then 00337 allocate(tvar2d_sin(sl%nx+1,sl%ny+1)) 00338 allocate(tvar2d_cos(sl%nx+1,sl%ny+1)) 00339 tvar2d_sin = 0.d0 00340 tvar2d_cos = 0.d0 00341 initialisedtvarsin = .true. 00342 endif 00343 00344 tvar2d_sin = tvar2d_sin + mult*sin(tvar2d) 00345 tvar2d_cos = tvar2d_cos + mult*cos(tvar2d) 00346 where (abs(tvar2d_sin)<numeps .and. abs(tvar2d_cos)<numeps) 00347 tvar2d_cos = numeps*sign(tvar2d_cos,1.d0) 00348 endwhere 00349 meansparslocal(i)%mean2d = atan2(tvar2d_sin,tvar2d_cos) 00350 !tvar2d_sin = nint(meansparsglobal(i)%mean2d) / 1d1 + nint(mult*sin(tvar2d)*1e6) 00351 !tvar2d_cos = mod(meansparsglobal(i)%mean2d,1.d0) * 1d7 + nint(mult*cos(tvar2d)*1e6) 00352 !meansparslocal(i)%mean2d = tvar2d_sin*1e1 + tvar2d_cos/1e7 00353 else 00354 meansparslocal(i)%mean2d = meansparslocal(i)%mean2d + mult*tvar2d 00355 endif 00356 meansparslocal(i)%variancecrossterm2d = & 00357 meansparslocal(i)%variancecrossterm2d/oldmean2d*meansparslocal(i)%mean2d + & 00358 mult*2.d0*tvar2d*meansparslocal(i)%mean2d 00359 meansparslocal(i)%variancesquareterm2d = & 00360 meansparslocal(i)%variancesquareterm2d+mult*(tvar2d)**2 00361 meansparslocal(i)%variance2d = & 00362 meansparslocal(i)%variancesquareterm2d-meansparslocal(i)%variancecrossterm2d + & 00363 meansparslocal(i)%mean2d**2 00364 meansparslocal(i)%max2d = max(meansparslocal(i)%max2d,tvar2d) 00365 meansparslocal(i)%min2d = min(meansparslocal(i)%min2d,tvar2d) 00366 case (3) 00367 allocate(oldmean3d(size(t%r3,1),size(t%r3,2),size(t%r3,3))) 00368 allocate(tvar3d(size(t%r3,1),size(t%r3,2),size(t%r3,3))) 00369 oldmean3d=meansparslocal(i)%mean3d 00370 ! bug in elsewere --> memory leak see commet Robert above 00371 where (oldmean3d<epsilon(0.d0) .and. oldmean3d>=0.d0) 00372 oldmean3d=epsilon(0.d0) 00373 endwhere 00374 where (oldmean3d>-epsilon(0.d0) .and. oldmean3d<0.d0) 00375 oldmean3d=-epsilon(0.d0) 00376 endwhere 00377 call gridrotate(par, sl,t,tvar3d) 00378 meansparslocal(i)%mean3d = meansparslocal(i)%mean3d + mult*tvar3d 00379 meansparslocal(i)%variancecrossterm3d = & 00380 meansparslocal(i)%variancecrossterm3d/oldmean3d*meansparslocal(i)%mean3d + & 00381 mult*2.d0*tvar3d*meansparslocal(i)%mean3d 00382 meansparslocal(i)%variancesquareterm3d = & 00383 meansparslocal(i)%variancesquareterm3d+mult*(tvar3d)**2 00384 meansparslocal(i)%variance3d = & 00385 meansparslocal(i)%variancesquareterm3d-meansparslocal(i)%variancecrossterm3d + & 00386 meansparslocal(i)%mean3d**2 00387 meansparslocal(i)%max3d = max(meansparslocal(i)%max3d,tvar3d) 00388 meansparslocal(i)%min3d = min(meansparslocal(i)%min3d,tvar3d) 00389 deallocate(oldmean3d,tvar3d) 00390 case (4) 00391 allocate(oldmean4d(size(t%r3,1),size(t%r3,2),size(t%r3,3),size(t%r4,4))) 00392 allocate(tvar4d(size(t%r3,1),size(t%r3,2),size(t%r3,3),size(t%r4,4))) 00393 oldmean4d=meansparslocal(i)%mean4d 00394 where (oldmean4d<epsilon(0.d0) .and. oldmean4d>=0.d0) 00395 oldmean4d=epsilon(0.d0) 00396 endwhere 00397 where (oldmean4d>-epsilon(0.d0) .and. oldmean4d<0.d0) 00398 oldmean4d=-epsilon(0.d0) 00399 endwhere 00400 call gridrotate(t,tvar4d) ! wwvv-todo 00401 meansparslocal(i)%mean4d = meansparslocal(i)%mean4d + mult*tvar4d 00402 meansparslocal(i)%variancecrossterm4d = & 00403 meansparslocal(i)%variancecrossterm4d/oldmean4d*meansparslocal(i)%mean4d + & 00404 mult*2.d0*tvar4d*meansparslocal(i)%mean4d 00405 meansparslocal(i)%variancesquareterm4d = & 00406 meansparslocal(i)%variancesquareterm4d+mult*(tvar4d)**2 00407 meansparslocal(i)%variance4d = & 00408 meansparslocal(i)%variancesquareterm4d-meansparslocal(i)%variancecrossterm4d + & 00409 meansparslocal(i)%mean4d**2 00410 meansparslocal(i)%max4d = max(meansparslocal(i)%max4d,tvar4d) 00411 meansparslocal(i)%min4d = min(meansparslocal(i)%min4d,tvar4d) 00412 deallocate(oldmean4d,tvar4d) 00413 end select 00414 enddo ! par%nmeanvar 00415 00416 end subroutine makeaverage 00417 00418 ! clear averages for new averaging period 00419 subroutine clearaverage(par) 00420 00421 use params 00422 00423 implicit none 00424 00425 integer :: i 00426 type(parameters), intent(in) :: par 00427 00428 do i=1,par%nmeanvar 00429 select case (meansparsglobal(i)%rank) 00430 case (2) 00431 meansparslocal(i)%mean2d = 0.d0 00432 meansparslocal(i)%variance2d = 0.d0 00433 meansparslocal(i)%variancecrossterm2d = 0.d0 00434 meansparslocal(i)%variancesquareterm2d = 0.d0 00435 meansparslocal(i)%min2d = huge(0.d0) 00436 meansparslocal(i)%max2d = -1.d0*huge(0.d0) 00437 case (3) 00438 meansparslocal(i)%mean3d = 0.d0 00439 meansparslocal(i)%variance3d = 0.d0 00440 meansparslocal(i)%variancecrossterm3d = 0.d0 00441 meansparslocal(i)%variancesquareterm3d = 0.d0 00442 meansparslocal(i)%min3d = huge(0.d0) 00443 meansparslocal(i)%max3d = -1.d0*huge(0.d0) 00444 case (4) 00445 meansparslocal(i)%mean4d = 0.d0 00446 meansparslocal(i)%variance4d = 0.d0 00447 meansparslocal(i)%variancecrossterm4d = 0.d0 00448 meansparslocal(i)%variancesquareterm4d = 0.d0 00449 meansparslocal(i)%min4d = huge(0.d0) 00450 meansparslocal(i)%max4d = -1.d0*huge(0.d0) 00451 end select 00452 enddo 00453 end subroutine clearaverage 00454 00455 00456 subroutine makecrossvector(s,sl,par,crossvararray,nvar,varindexvec,mg,cstype) 00457 use params, only: parameters 00458 use spaceparams 00459 use mnemmodule 00460 00461 IMPLICIT NONE 00462 00463 type(spacepars), intent(in) :: s,sl 00464 type(parameters) :: par 00465 real*8,dimension(:,:) :: crossvararray 00466 integer,intent(in) :: mg,cstype,nvar 00467 integer,dimension(nvar),intent(in):: varindexvec 00468 type(arraytype) :: t 00469 00470 integer :: i 00471 #ifdef USEMPI 00472 logical :: toall = .true. 00473 #endif 00474 00475 crossvararray=-999.d0 00476 ! wwvv to avoid warning about unused sl and par 00477 if (.false.) then 00478 print *,sl%nx,par%swave 00479 endif 00480 do i=1,nvar 00481 #ifdef USEMPI 00482 call space_collect_index(s,sl,par,varindexvec(i)) 00483 #endif 00484 if (xmaster) then 00485 call indextos(s,varindexvec(i),t) 00486 select case(t%type) 00487 case('r') 00488 if(cstype==0) then 00489 crossvararray(i,:)=t%r2(:,mg) 00490 else 00491 crossvararray(i,:)=t%r2(mg,:) 00492 endif 00493 case('i') 00494 if(cstype==0) then 00495 crossvararray(i,:)=t%i2(:,mg) 00496 else 00497 crossvararray(i,:)=t%i2(mg,:) 00498 endif 00499 end select 00500 endif 00501 enddo 00502 #ifdef USEMPI 00503 call xmpi_bcast(crossvararray,toall) ! wwvv todo: it seems that crosvararry is only needed at xomaster 00504 #endif 00505 ! wwvv to avoid warning about unused sl: 00506 if (sl%nx .eq. -1) return 00507 end subroutine makecrossvector 00508 00509 #ifdef USEMPI 00510 subroutine means_collect(sl,a,b) 00511 ! 00512 ! collect mean variables, output in a 00513 ! 00514 use spaceparams 00515 00516 implicit none 00517 00518 type(spacepars),intent(in) :: sl 00519 type(meanspars),intent(inout) :: a 00520 type(meanspars),intent(in) :: b 00521 00522 select case (b%rank) 00523 case (2) 00524 call space_collect(sl,a%mean2d,b%mean2d) 00525 call space_collect(sl,a%variance2d,b%variance2d) 00526 call space_collect(sl,a%max2d,b%max2d) 00527 call space_collect(sl,a%min2d,b%min2d) 00528 ! for H and urms : 00529 call space_collect(sl,a%variancesquareterm2d,b%variancesquareterm2d) 00530 case (3) 00531 call space_collect(sl,a%mean3d,b%mean3d) 00532 call space_collect(sl,a%variance3d,b%variance3d) 00533 call space_collect(sl,a%max3d,b%max3d) 00534 call space_collect(sl,a%min3d,b%min3d) 00535 case (4) 00536 call space_collect(sl,a%mean4d,b%mean4d) 00537 call space_collect(sl,a%variance4d,b%variance4d) 00538 call space_collect(sl,a%max4d,b%max4d) 00539 call space_collect(sl,a%min4d,b%min4d) 00540 end select 00541 00542 end subroutine means_collect 00543 #endif 00544 00545 00546 end module means_module