XBeach
|
00001 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00002 ! wwvv-todo many changes wrt trunk 00003 !!!!!!!!!!!!!!!!!!!!!!! MODULE POSTPROCESS !!!!!!!!!!!!!!!!!!!!!!! 00004 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00005 00006 ! This module contains transformation routines from xbeach internal data structures to output data structures. 00007 00008 module postprocessmod 00009 implicit none 00010 save 00011 !interface postprocessvar 00012 ! module procedure postprocessvar_r2 00013 !end interface postprocessvar 00014 interface gridrotate 00015 ! rotate grids, given s, t and outputs x 00016 module procedure gridrotate_r0 00017 module procedure gridrotate_r1 00018 module procedure gridrotate_r2 00019 module procedure gridrotate_r3 00020 module procedure gridrotate_r4 00021 module procedure gridrotate_i0 00022 module procedure gridrotate_i1 00023 module procedure gridrotate_i2 00024 module procedure gridrotate_i3 00025 module procedure gridrotate_i4 00026 end interface gridrotate 00027 00028 contains 00029 subroutine postprocessvar_r2(wetz, varname, r2, dFill) 00030 use spaceparams 00031 use mnemmodule 00032 00033 integer, dimension(:,:), intent(in) :: wetz 00034 character(maxnamelen),intent(in) :: varname 00035 real*8, intent(in) :: dFill 00036 real*8, dimension(:,:), intent(inout) :: r2 00037 00038 select case(varname) 00039 case(mnem_zs,mnem_u,mnem_v,mnem_H,mnem_thetamean,mnem_k,mnem_sigm,mnem_c,mnem_cg,mnem_hh,mnem_E,mnem_R, & 00040 mnem_urms,mnem_D,mnem_ue,mnem_ve,mnem_Sk,mnem_As) 00041 where (wetz==0) 00042 r2 = dFill 00043 elsewhere 00044 r2 = r2 00045 endwhere 00046 case default 00047 ! do nothing 00048 end select 00049 00050 end subroutine postprocessvar_r2 00051 00052 subroutine snappointstogrid(par, s, xpoints, ypoints,scrprintinp,dmin) 00053 ! Lookup the nearest neighbour grid coordinates of the output points specified in params file 00054 ! Convert world coordinates of points to nearest (lsm) grid point 00055 use spaceparams 00056 use params 00057 use logging_module 00058 00059 type(spacepars), intent(in) :: s 00060 type(parameters), intent(in) :: par 00061 00062 integer*4,dimension(:),intent(inout) :: xpoints ! model x-coordinate of output points 00063 integer*4,dimension(:),intent(inout) :: ypoints ! model y-coordinate of output points 00064 00065 real*8,dimension(s%nx+1,s%ny+1) :: mindist 00066 integer,dimension(2) :: minlocation 00067 integer :: i 00068 logical,optional,intent(in) :: scrprintinp 00069 logical :: scrprint 00070 real*8,dimension(:),allocatable :: mindistr 00071 real*8,dimension(:),optional,intent(out) :: dmin 00072 character(1000) :: txt 00073 ! 00074 if (present(scrprintinp)) then 00075 scrprint=scrprintinp 00076 else 00077 scrprint=.true. 00078 endif 00079 ! 00080 allocate(mindistr(par%npoints+par%nrugauge)) 00081 ! Let's hope that the s%xz and s%yz are already available.... 00082 ! Compute the minimum distances for each point 00083 if (par%npoints + par%nrugauge > 0) then 00084 do i=1,(par%npoints+par%nrugauge) 00085 mindist=sqrt((par%xpointsw(i)-s%xz)**2+(par%ypointsw(i)-s%yz)**2) 00086 ! look up the location of the found minimum 00087 minlocation=minloc(mindist) 00088 ! minimum distance 00089 mindistr(i) = mindist(minlocation(1),minlocation(2)) 00090 ! The y coordinate is always the same 00091 ypoints(i)=minlocation(2) 00092 ! For rugauges the xpoint is always 1 00093 if (par%pointtypes(i) == 1) then 00094 xpoints(i)=1 00095 if (scrprint) call writelog('ls','(a,i0)','Runup gauge at grid line iy=',ypoints(i)) 00096 else 00097 xpoints(i)=minlocation(1) 00098 if(scrprint) then ! wwvv-todo 00099 write(txt,"('Distance output point',i3.3,'(',f0.2,',',f0.2,') to gridpoint(',i0,',',i0'): ',f0.2,' m')") & 00100 i,par%xpointsw(i),par%ypointsw(i),xpoints(i),ypoints(i),mindistr(i) 00101 call writelog('ls','(a)',txt) !,xomaster) 00102 endif 00103 endif 00104 end do 00105 end if 00106 if (present(dmin)) then 00107 dmin = mindistr 00108 endif 00109 deallocate(mindistr) 00110 end subroutine snappointstogrid 00111 00112 00113 subroutine gridrotate_r0(t, x) 00114 use params 00115 use spaceparams 00116 use mnemmodule 00117 00118 type(arraytype), intent(in) :: t 00119 real*8 :: x 00120 x = t%r0 00121 end subroutine gridrotate_r0 00122 00123 subroutine gridrotate_r1(par, s, t, x) 00124 use params 00125 use spaceparams 00126 use mnemmodule 00127 00128 type(parameters), intent(in) :: par 00129 type(spacepars), intent(in) :: s 00130 type(arraytype), intent(in) :: t 00131 real*8, dimension(:) :: x 00132 real*8, parameter :: pi = 4*atan(1.0d0) 00133 00134 if (par%rotate .eq. 1 .and. s%nx .ne. -1) then 00135 select case(t%name) 00136 case(mnem_theta) 00137 x=270-(t%r1*(180/pi)) 00138 case(mnem_theta0) 00139 x=270-(t%r1*(180/pi)) 00140 case default 00141 x=t%r1 00142 end select 00143 else 00144 x=t%r1 00145 endif 00146 end subroutine gridrotate_r1 00147 00148 00149 subroutine gridrotate_r2(par, s, t, x , row, col ) 00150 ! 00151 ! sl: if present then this subroutine needs to be called by everyone 00152 ! in xmpi_ocomm, the extra needed s% variables needed will be collected 00153 ! x will not assigned to 00154 ! 00155 ! if row is present, and row>0, the rotation will be done only for 00156 ! s%some_arry(row,col) and assigned to x(1,1) 00157 00158 use params 00159 use spaceparams 00160 use mnemmodule 00161 00162 type(parameters), intent(in) :: par 00163 type(spacepars), intent(inout) :: s 00164 type(arraytype), intent(in) :: t 00165 real*8, dimension(:,:), intent(out) :: x 00166 integer, intent(in), optional :: row, col 00167 00168 real*8, parameter :: pi = 4*atan(1.0d0) 00169 real*8, dimension(size(s%alfaz,1), size(s%alfaz,2)) :: Sutot, Svtot 00170 00171 integer is, ie, js, je 00172 integer isx, iex, jsx, jex 00173 00174 Sutot = 0; 00175 Svtot = 0; 00176 00177 is = 1 00178 ie = size(x,1) 00179 js = 1 00180 je = size(x,2) 00181 isx = is 00182 iex = ie 00183 jsx = js 00184 jex = je 00185 00186 if (present(row)) then 00187 if (row .gt. 0 ) then 00188 is = row 00189 ie = row 00190 js = col 00191 je = col 00192 isx = 1 00193 iex = 1 00194 jsx = 1 00195 jex = 1 00196 endif 00197 endif 00198 00199 if (par%rotate .eq. 1) then 00200 select case(t%name) 00201 case(mnem_thetamean) 00202 x(isx:iex,jsx:jex) = 270-((t%r2(is:ie,js:je))*(180/pi)) 00203 case(mnem_Fx) 00204 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00205 & s%Fy(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00206 case(mnem_Fy) 00207 x(isx:iex,jsx:jex) = s%Fx(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00208 & t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00209 case(mnem_u) 00210 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00211 & s%v (is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00212 case(mnem_gwu) 00213 x(isx:iex,jsx:jex) = t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00214 & s%gwv(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00215 case(mnem_v) 00216 x(isx:iex,jsx:jex) = s%u (is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00217 & t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00218 case(mnem_gwv) 00219 x(isx:iex,jsx:jex) = s%gwu(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00220 & t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00221 case(mnem_ue) 00222 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00223 & s%ve(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00224 case(mnem_ve) 00225 x(isx:iex,jsx:jex) = s%ue(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00226 & t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00227 case(mnem_ui) 00228 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00229 & s%vi(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00230 case(mnem_vi) 00231 x(isx:iex,jsx:jex) = s%ui(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00232 & t%r2(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00233 case(mnem_umean) 00234 x(isx:iex,jsx:jex) = t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00235 & s%vmean(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00236 case(mnem_vmean) 00237 x(isx:iex,jsx:jex) = s%umean(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00238 & t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00239 case(mnem_uwf) 00240 x(isx:iex,jsx:jex) = t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je))- & 00241 & s%vwf(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00242 case(mnem_vwf) 00243 x(isx:iex,jsx:jex) = s%uwf(is:ie,js:je)*sin(s%alfaz(is:ie,js:je))+ & 00244 & t%r2 (is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00245 case(mnem_Sutot) 00246 ! Jaap interpolate transports to water level points before rotating to real world 00247 Sutot = sum(s%Subg,dim=3)+sum(s%Susg,dim=3) 00248 Svtot = sum(s%Svbg,dim=3)+sum(s%Svsg,dim=3) 00249 Sutot(2:s%nx,:)=0.5d0*(Sutot(1:s%nx-1,:)+Sutot(2:s%nx,:)) 00250 ! Dano: take care of sf1d! 00251 if (s%ny>0) then 00252 Svtot(:,2:s%ny)=0.5d0*(Svtot(:,1:s%ny-1)+Svtot(:,2:s%ny)) 00253 endif 00254 x(isx:iex,jsx:jex) = Sutot(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) - & 00255 & Svtot(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) 00256 case(mnem_Svtot) 00257 Sutot = sum(s%Subg,dim=3)+sum(s%Susg,dim=3) 00258 Svtot = sum(s%Svbg,dim=3)+sum(s%Svsg,dim=3) 00259 Sutot(2:s%nx,:)=0.5d0*(Sutot(1:s%nx-1,:)+Sutot(2:s%nx,:)) 00260 Svtot(:,2:s%ny)=0.5d0*(Svtot(:,1:s%ny-1)+Svtot(:,2:s%ny)) 00261 x(isx:iex,jsx:jex) = Sutot(is:ie,js:je)*sin(s%alfaz(is:ie,js:je)) + & 00262 & Svtot(is:ie,js:je)*cos(s%alfaz(is:ie,js:je)) 00263 case(mnem_cctot) 00264 x(isx:iex,jsx:jex) = sum(s%ccg(is:ie,js:je,:),dim=3) 00265 case default 00266 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je) 00267 end select 00268 else 00269 select case(t%name) 00270 case(mnem_Sutot) 00271 Sutot = sum(s%Subg,dim=3)+sum(s%Susg,dim=3) 00272 Sutot(2:s%nx,:)=0.5d0*(Sutot(1:s%nx-1,:)+Sutot(2:s%nx,:)) 00273 case(mnem_Svtot) 00274 Sutot = sum(s%Subg,dim=3)+sum(s%Susg,dim=3) 00275 Svtot(:,2:s%ny)=0.5d0*(Svtot(:,1:s%ny-1)+Svtot(:,2:s%ny)) 00276 case default 00277 x(isx:iex,jsx:jex) = t%r2(is:ie,js:je) 00278 end select 00279 endif 00280 00281 end subroutine gridrotate_r2 00282 00283 subroutine gridrotate_r3(par, s, t, x, row, col ) 00284 ! if row is present, and row>0, the rotation will be done only for 00285 ! s%some_arry(row,col,:) and assigned to x(1,1,:) 00286 00287 use params 00288 use spaceparams 00289 use mnemmodule 00290 use logging_module 00291 00292 implicit none 00293 00294 type(parameters), intent(in) :: par 00295 type(spacepars), intent(in) :: s 00296 type(arraytype), intent(in) :: t 00297 real*8, dimension(:,:,:) :: x 00298 integer, intent(in), optional :: row, col 00299 00300 real*8, dimension(size(s%alfaz,1), size(s%alfaz,2), size(t%r3,3)) :: alfazr3,Susg,Svsg,Subg,Svbg 00301 real*8, parameter :: pi = 4*atan(1.0d0) 00302 integer :: i 00303 00304 integer is, ie, js, je 00305 integer isx, iex, jsx, jex 00306 00307 if (present(row)) then 00308 if (row .gt. 0 ) then 00309 is = row 00310 ie = row 00311 js = col 00312 je = col 00313 isx = 1 00314 iex = 1 00315 jsx = 1 00316 jex = 1 00317 endif 00318 else 00319 is = 1 00320 ie = size(x,1) 00321 js = 1 00322 je = size(x,2) 00323 isx = is 00324 iex = ie 00325 jsx = js 00326 jex = je 00327 endif 00328 00329 ! Jaap: initialize local transport variables (used for interpolation to water level points) 00330 Susg = 0; 00331 Svsg = 0; 00332 Subg = 0; 00333 Svbg = 0; 00334 00335 ! fill variable alfazr3. We presume first 2 dimension are related to nx+1, ny+1 respectively 00336 ! Better double check 00337 if (size(s%alfaz,1) .ne. size(t%r3,1)) call writelog('ls', '', & 00338 & 'Assertion error, s%alfaz and t%r3 do not align on 1st dimension', size(s%alfaz,1), size(t%r3,1)) 00339 if (size(s%alfaz,2) .ne. size(t%r3,2)) call writelog('ls', '', & 00340 & 'Assertion error, s%alfaz and t%r3 do not align on 2nd dimension', size(s%alfaz,2), size(t%r3,2)) 00341 ! This should be something like: 00342 ! alfazr3 = (/(s%alfaz, i=1,size(t%r3,3)) /) 00343 do i=1,size(t%r3,3) 00344 alfazr3(:,:,i) = s%alfaz 00345 end do 00346 00347 if (par%rotate .eq. 1) then 00348 select case(t%name) 00349 case(mnem_cgx) 00350 x(isx:iex,jsx:jex,:) = t%r3 (is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:))- & 00351 & s%cgy(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:)) 00352 case(mnem_cgy) 00353 x(isx:iex,jsx:jex,:) = s%cgx(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:))+ & 00354 & t%r3 (is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:)) 00355 case(mnem_cx) 00356 x(isx:iex,jsx:jex,:) = t%r3(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:))- & 00357 & s%cy(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:)) 00358 case(mnem_cy) 00359 x(isx:iex,jsx:jex,:) = s%cx(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:))+ & 00360 & t%r3(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:)) 00361 case(mnem_thet) 00362 x(isx:iex,jsx:jex,:) = 270-((s%thet(is:ie,js:je,:)+alfazr3(is:ie,js:je,:))*(180/pi)) 00363 case(mnem_Susg) 00364 ! Jaap: interpolate transports to water level points before rotating to world coordinates 00365 Susg = t%r3 00366 Svsg = s%Svsg 00367 Susg(2:s%nx,:,:)=0.5d0*(t%r3(1:s%nx-1,:,:)+t%r3(2:s%nx,:,:)) 00368 if (s%ny>0) then 00369 Svsg(:,2:s%ny,:)=0.5d0*(s%Svsg(:,1:s%ny-1,:)+s%Svsg(:,2:s%ny,:)) 00370 endif 00371 x(isx:iex,jsx:jex,:) = Susg(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:))- & 00372 & Svsg(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:)) 00373 case(mnem_Svsg) 00374 Susg = s%Susg 00375 Svsg = t%r3 00376 Susg(2:s%nx,:,:)=0.5d0*(s%Susg(1:s%nx-1,:,:)+s%Susg(2:s%nx,:,:)) 00377 if (s%ny>0) then 00378 Svsg(:,2:s%ny,:)=0.5d0*(t%r3(:,1:s%ny-1,:)+t%r3(:,2:s%ny,:)) 00379 endif 00380 x(isx:iex,jsx:jex,:) = Susg(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:))+ & 00381 & Svsg(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:)) 00382 case(mnem_Subg) 00383 Subg = t%r3 00384 Svbg = s%Svbg 00385 Subg(2:s%nx,:,:)=0.5d0*(t%r3(1:s%nx-1,:,:)+t%r3(2:s%nx,:,:)) 00386 if (s%ny>0) then 00387 Svbg(:,2:s%ny,:)=0.5d0*(s%Svbg(:,1:s%ny-1,:)+s%Svbg(:,2:s%ny,:)) 00388 endif 00389 x(isx:iex,jsx:jex,:) = Subg(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:))- & 00390 & Svbg(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:)) 00391 case(mnem_Svbg) 00392 Subg = s%Subg 00393 Svbg = t%r3 00394 Subg(2:s%nx,:,:)=0.5d0*(s%Subg(1:s%nx-1,:,:)+s%Subg(2:s%nx,:,:)) 00395 if (s%ny>0) then 00396 Svbg(:,2:s%ny,:)=0.5d0*(t%r3(:,1:s%ny-1,:)+t%r3(:,2:s%ny,:)) 00397 endif 00398 x(isx:iex,jsx:jex,:) = Subg(is:ie,js:je,:)*sin(alfazr3(is:ie,js:je,:))+ & 00399 & Svbg(is:ie,js:je,:)*cos(alfazr3(is:ie,js:je,:)) 00400 case default 00401 x(isx:iex,jsx:jex,:) = t%r3(is:ie,js:je,:) 00402 end select 00403 else 00404 00405 x(isx:iex,jsx:jex,:) = t%r3(is:ie,js:je,:) 00406 endif 00407 end subroutine gridrotate_r3 00408 00409 subroutine gridrotate_r4(t, x, row, col) 00410 use params 00411 use spaceparams 00412 use mnemmodule 00413 00414 type(arraytype), intent(in) :: t 00415 real*8, dimension(:,:,:,:) :: x 00416 integer, optional, intent(in) :: row,col 00417 integer is, ie, js, je 00418 integer isx, iex, jsx, jex 00419 00420 if (present(row)) then 00421 if (row .gt. 0 ) then 00422 is = row 00423 ie = row 00424 js = col 00425 je = col 00426 isx = 1 00427 iex = 1 00428 jsx = 1 00429 jex = 1 00430 endif 00431 else 00432 is = 1 00433 ie = size(x,1) 00434 js = 1 00435 je = size(x,2) 00436 isx = is 00437 iex = ie 00438 jsx = js 00439 jex = je 00440 endif 00441 00442 x(isx:iex,jsx:jex,:,:) = t%r4(is:ie,js:je,:,:) 00443 00444 end subroutine gridrotate_r4 00445 00446 00447 subroutine gridrotate_i0(t, x) 00448 use params 00449 use spaceparams 00450 use mnemmodule 00451 00452 type(arraytype), intent(in) :: t 00453 integer :: x 00454 x = t%i0 00455 end subroutine gridrotate_i0 00456 00457 subroutine gridrotate_i1(t, x) 00458 use params 00459 use spaceparams 00460 use mnemmodule 00461 00462 type(arraytype), intent(in) :: t 00463 integer, dimension(:),intent(out) :: x 00464 x = t%i1 00465 end subroutine gridrotate_i1 00466 00467 subroutine gridrotate_i2(t, x) 00468 use params 00469 use spaceparams 00470 use mnemmodule 00471 00472 type(arraytype), intent(in) :: t 00473 integer, dimension(:,:),intent(out) :: x 00474 x = t%i2 00475 end subroutine gridrotate_i2 00476 00477 subroutine gridrotate_i3(t, x) 00478 use params 00479 use spaceparams 00480 use mnemmodule 00481 00482 type(arraytype), intent(in) :: t 00483 integer, dimension(:,:,:),intent(out) :: x 00484 x = t%i3 00485 end subroutine gridrotate_i3 00486 00487 subroutine gridrotate_i4(t, x) 00488 use params 00489 use spaceparams 00490 use mnemmodule 00491 00492 type(arraytype), intent(in) :: t 00493 integer, dimension(:,:,:,:),intent(out) :: x 00494 x = t%i4 00495 end subroutine gridrotate_i4 00496 00497 00498 function runup(par, s) 00499 use spaceparams 00500 use params 00501 type(spacepars) :: s 00502 type(parameters) :: par 00503 integer, dimension(par%npoints + par%nrugauge) :: runup 00504 00505 ! This is not correct.... 00506 where (par%pointtypes == 1) 00507 runup = max(maxval(minloc(s%wetz(:,:)))-1,1) 00508 end where 00509 end function runup 00510 00511 function get_sister_mnem(mnem) result(sistermnem) 00512 use mnemmodule 00513 character(maxnamelen) :: sistermnem 00514 character(maxnamelen) :: mnem 00515 00516 select case (mnem) 00517 case(mnem_Fx) 00518 sistermnem = mnem_Fy 00519 case(mnem_Fy) 00520 sistermnem = mnem_Fx 00521 case(mnem_u) 00522 sistermnem = mnem_v 00523 case(mnem_gwu) 00524 sistermnem = mnem_gwv 00525 case(mnem_v) 00526 sistermnem = mnem_u 00527 case(mnem_gwv) 00528 sistermnem = mnem_gwu 00529 case(mnem_ue) 00530 sistermnem = mnem_ve 00531 case(mnem_ve) 00532 sistermnem = mnem_ue 00533 case(mnem_ui) 00534 sistermnem = mnem_vi 00535 case(mnem_vi) 00536 sistermnem = mnem_ui 00537 case(mnem_umean) 00538 sistermnem = mnem_vmean 00539 case(mnem_vmean) 00540 sistermnem = mnem_umean 00541 case(mnem_uwf) 00542 sistermnem = mnem_vwf 00543 case(mnem_vwf) 00544 sistermnem = mnem_uwf 00545 case(mnem_Sutot) 00546 sistermnem = mnem_Svtot 00547 case(mnem_Svtot) 00548 sistermnem = mnem_Sutot 00549 case(mnem_cctot) 00550 sistermnem = mnem_ccg 00551 case(mnem_cgx) 00552 sistermnem = mnem_cgy 00553 case(mnem_cgy) 00554 sistermnem = mnem_cgx 00555 case(mnem_cx) 00556 sistermnem = mnem_cy 00557 case(mnem_cy) 00558 sistermnem = mnem_cx 00559 case(mnem_Susg) 00560 sistermnem = mnem_Svsg 00561 case(mnem_Svsg) 00562 sistermnem = mnem_Susg 00563 case(mnem_Subg) 00564 sistermnem = mnem_Svbg 00565 case(mnem_Svbg) 00566 sistermnem = mnem_Subg 00567 case default 00568 sistermnem = 'none' 00569 end select 00570 end function get_sister_mnem 00571 end module postprocessmod