!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !$Author: joelevin $ !$Date: 2011-03-24 11:33:26 -0500 (Thu, 24 Mar 2011) $ !$Revision: 11724 $ !$HeadURL: https://svn.weru.ksu.edu/weru/weps1/trunk/weps.src/src/lib_erosion/daily_erodout.for $ subroutine daily_erodout (o_unit, o_E_unit) use p1werm_def use datetime_def use m1geo_def use gridmod IMPLICIT NONE ! +++ PURPOSE +++ ! To print output desired from standalone EROSION submodel ! +++ ARGUMENT DECLARATIONS +++ integer o_unit, o_E_unit integer i, j real asoilLossTot, asoilLossSusp, asoilLossPM10 real tt, lx, ly real topt,topss, top10, bott, botss, bot10, ritt, ritss, rit10 real lftt, lftss, lft10, tot, totbnd integer initflag save initflag integer yr, mo, da ! + + + GLOBAL COMMON BLOCKS + + + ! include 'p1werm.inc' include 'm1flag.inc' ! + + + LOCAL COMMON BLOCKS ! include 'erosion/e2erod.inc' include 'erosion/m2geo.inc' ! integer x, y ! ++++ ARGUMENT DEFINITIONS +++ ! +++ SUBROUTINES CALLED+++ ! plotout.for ! ++++ LOCAL VARIABLES +++ ! +++ END SPECIFICATIONS +++ ! Calculate Averages Crossing Borders ! top border asoilLossTot = 0.0 asoilLossSusp = 0.0 asoilLossPM10 = 0.0 j = jmax do 1 i = 1, imax-1 asoilLossTot = asoilLossTot + grid(i,j)%soilLossTot asoilLossSusp = asoilLossSusp + grid(i,j)%soilLossSusp asoilLossPM10 = asoilLossPM10 + grid(i,j)%soilLossPM10 1 continue ! calc. average at top border topt = asoilLossTot/(imax-1) topss = asoilLossSusp/(imax-1) top10 = asoilLossPM10/(imax-1) ! bottom border asoilLossTot = 0.0 asoilLossSusp = 0.0 asoilLossPM10 = 0.0 j = 0 do 2 i = 1, imax-1 asoilLossTot = asoilLossTot + grid(i,j)%soilLossTot asoilLossSusp = asoilLossSusp + grid(i,j)%soilLossSusp asoilLossPM10 = asoilLossPM10 + grid(i,j)%soilLossPM10 2 continue ! calc. average at bottom border bott = asoilLossTot/(imax-1) botss = asoilLossSusp/(imax-1) bot10 = asoilLossPM10/(imax-1) ! right border asoilLossTot = 0.0 asoilLossSusp = 0.0 asoilLossPM10 = 0.0 i = imax do 3 j = 1, jmax-1 asoilLossTot = asoilLossTot + grid(i,j)%soilLossTot asoilLossSusp = asoilLossSusp + grid(i,j)%soilLossSusp asoilLossPM10 = asoilLossPM10 + grid(i,j)%soilLossPM10 3 continue ! calc. average at right border ritt = asoilLossTot/(jmax-1) ritss = asoilLossSusp/(jmax-1) rit10 = asoilLossPM10/(jmax-1) ! ! left border asoilLossTot = 0.0 asoilLossSusp = 0.0 asoilLossPM10 = 0.0 i = 0 do 4 j = 1, jmax-1 asoilLossTot = asoilLossTot + grid(i,j)%soilLossTot asoilLossSusp = asoilLossSusp + grid(i,j)%soilLossSusp asoilLossPM10 = asoilLossPM10 + grid(i,j)%soilLossPM10 4 continue ! calc. average at left border lftt = asoilLossTot/(jmax-1) lftss = asoilLossSusp/(jmax-1) lft10 = asoilLossPM10/(jmax-1) ! calculate averages of inner grid points asoilLossTot = 0.0 asoilLossSusp = 0.0 asoilLossPM10 = 0.0 do 5 j=1,jmax-1 do 5 i= 1, imax-1 asoilLossTot= asoilLossTot + grid(i,j)%soilLossTot asoilLossSusp = asoilLossSusp + grid(i,j)%soilLossSusp asoilLossPM10 = asoilLossPM10 + grid(i,j)%soilLossPM10 5 continue tt = (imax-1)*(jmax-1) asoilLossTot = asoilLossTot/tt asoilLossSusp = asoilLossSusp/tt asoilLossPM10 = asoilLossPM10/tt ! calculate comparision of boundary and interior losses lx = amxsim(1,2) - amxsim(1,1) ly = amxsim(2,2) - amxsim(2,1) tot = asoilLossTot*lx*ly totbnd = (topt + bott + topss + botss)*lx + & & (ritt + lftt + ritss + lftss)*ly if (btest(am0efl,1)) then if (initflag .eq. 0) then !write header to files write (o_unit,*) 'Grid cell output from daily_erodout.for' write (o_unit,*) ! Print date of Run write (o_unit,*) 'Date of run: ', datetimestr write (o_unit,*) write (o_unit,fmt="(1x,a,5f10.2)") "", & & amasim,((amxsim(x,y),x=1,2),y=1,2) write (o_unit,*) write (o_unit,*) & & 'Total grid size: (', imax+1,',', jmax+1, ') ',& & 'Inner grid size: (', imax-1,',', jmax-1, ')' write (o_unit,*) initflag = initflag + 1 endif call caldatw (da, mo, yr) !get day, month and year write (UNIT=o_unit,FMT="(' ',i5,i3,i3,' ')",ADVANCE="YES") & & yr, mo, da write (o_unit,*) write (o_unit,*) & & ' top(i=1,imax-1,j=jmax) ', & & 'bottom(i=1,imax-1,j=0) ', & & 'right(i=imax,j=1,jmax-1) ', & & 'left(i=0, j=1,jmax-1) ' write (o_unit,10) (grid(i,jmax)%soilLossTot+grid(i,jmax)%soilLossSusp, i = 1, imax-1) write (o_unit,10) (grid(i,0)%soilLossTot+grid(i,0)%soilLossSusp, i = 1, imax-1) write (o_unit,10) (grid(imax,j)%soilLossTot+grid(imax,j)%soilLossSusp, j = 1, jmax-1) write (o_unit,10) (grid(0,j)%soilLossTot+grid(0,j)%soilLossSusp, j = 1, jmax-1) write (o_unit,*) write (o_unit,7) write (o_unit,*) & & ' top(i=1,imax-1,j=jmax) ', & & 'bottom(i=1,imax-1,j=0) ', & & 'right(i=imax,j=1,jmax-1) ', & & 'left(i=0, j=1,jmax-1) ' write (o_unit,10) (grid(i,jmax)%soilLossTot, i = 1, imax-1) write (o_unit,10) (grid(i,0)%soilLossTot, i = 1, imax-1) write (o_unit,10) (grid(imax,j)%soilLossTot, j = 1, jmax-1) write (o_unit,10) (grid(0,j)%soilLossTot, j = 1, jmax-1) write (o_unit,*) write (o_unit,8) write (o_unit,*) & & ' top(i=1,imax-1,j=jmax) ', & & 'bottom(i=1,imax-1,j=0) ', & & 'right(i=imax,j=1,jmax-1) ', & & 'left(i=0, j=1,jmax-1) ' write (o_unit,10) (grid(i,jmax)%soilLossSusp, i = 1, imax-1) write (o_unit,10) (grid(i,0)%soilLossSusp, i = 1, imax-1) write (o_unit,10) (grid(jmax,j)%soilLossSusp, j = 1, jmax-1) write (o_unit,10) (grid(0,j)%soilLossSusp, j = 1, jmax-1) write (o_unit,*) write (o_unit,9) write (o_unit,*) & & ' top(i=1,imax-1,j=jmax) ', & & 'bottom(i=1,imax-1,j=0) ', & & 'right(i=imax,j=1,jmax-1) ', & & 'left(i=0,j=1,jmax-1) ' write (o_unit,11) (grid(i,jmax)%soilLossPM10, i = 1, imax-1) write (o_unit,11) (grid(i,0)%soilLossPM10, i = 1, imax-1) write (o_unit,11) (grid(imax,j)%soilLossPM10, j = 1, jmax-1) write (o_unit,11) (grid(0,j)%soilLossPM10, j = 1, jmax-1) write (o_unit,*) write (o_unit,fmt="(' | ',3('|',a))") & & 'Total Soil Loss', 'soil loss', '(kg/m^2)' do 19 j = jmax-1, 1, -1 write (o_unit,10) (grid(i,j)%soilLossTot, i = 1, imax-1) 19 continue write (o_unit,fmt="(' ')") write (o_unit,*) write (o_unit,fmt="(' | ',3('|',a))") & & 'Saltation/Creep Soil Loss', 'salt/creep soil loss', '(kg/m^2)' do 29 j = jmax-1, 1, -1 write (o_unit,10) (grid(i,j)%soilLossTot-grid(i,j)%soilLossSusp, i = 1, imax-1) 29 continue write (o_unit,fmt="(' ')") write (o_unit,*) write (o_unit,fmt="(' | ',3('|',a))") & & 'Suspension Soil Loss', 'suspension soil loss', '(kg/m^2)' do 39 j = jmax-1, 1, -1 write (o_unit,10) (grid(i,j)%soilLossSusp, i = 1, imax-1) 39 continue write (o_unit,fmt="(' ')") write (o_unit,*) write (o_unit,fmt="(' | ',3('|',a))") & & 'PM10 Soil Loss', 'PM10 soil loss', '(kg/m^2)' do 49 j = jmax-1, 1, -1 write (o_unit,11) (grid(i,j)%soilLossPM10, i = 1, imax-1) 49 continue write (o_unit,fmt="(' ')") write (o_unit,*) write (o_unit,*) '**Averages - Field' write (o_unit,*) ' Total salt/creep susp PM10 ' write (o_unit,*) ' soilLossTot soilLossSusp soilLossPM10' write (o_unit,*) ' -----------------kg/m^2--------------------' write (o_unit,15) asoilLossTot, asoilLossTot-asoilLossSusp, asoilLossSusp, asoilLossPM10 write (o_unit,*) write (o_unit,*) '**Averages - Crossing Boundaries ' write (o_unit,*) 'Location Total Salt/Creep Susp PM10' write (o_unit,*) '--------------------kg/m----------------------' write (o_unit,21) topt+topss, topt, topss, top10 write (o_unit,22) bott+botss, bott, botss, bot10 write (o_unit,23) ritt+ritss, ritt, ritss, rit10 write (o_unit,24) lftt+lftss, lftt, lftss, lft10 write (o_unit,*) write (o_unit,*) ' Comparision of interior & boundary loss' write (o_unit,*) ' interior boundary int/bnd ratio' if( totbnd.gt.1.0e-9 ) then write (o_unit,16) tot, totbnd, tot/totbnd else !Boundary loss near or equal to zero write (o_unit,16) tot, totbnd, 1.0e-9 end if !^^^tmp out ! write (o_unit,*) 'lx=', lx, 'ly=', ly,'tot',tot !^^^end tmpout ! additional output statements for easy shell script parsing write (o_unit,*) ! write deposition as positive numbers write (o_unit,17) asoilLossTot, asoilLossTot-asoilLossSusp, asoilLossSusp, asoilLossPM10 17 format (' repeat of total, salt/creep, susp, PM10:', 3f12.4,f12.6) ! output formats ! 6 format (1x,' Passing Border Grid Cells - Total soilLossTot+soilLossSusp(kg/m)') 7 format (1x,' Passing Border Grid Cells - Salt/Creep soilLossTot(kg/m)') 8 format (1x,' Passing Border Grid Cells - Suspension soilLossSusp(kg/m)') 9 format (1x,' Passing Border Grid Cells - PM10 soilLossPM10(kg/m)') ! 50 format (1x,' Leaving Field Grid Cells - Total soilLossTot(kg/m^2)') ! 60 format (1x,' Leaving Field Grid Cells - Salt/Creep soilLossTot-soilLossSusp(kg/m& ! &^2)') ! 70 format (1x,' Leaving Field Grid Cells - Suspension soilLossSusp(kg/m^2) & ! &') ! 80 format (1x,' Leaving Field Grid Cells - PM10 soilLossPM10(kg/m^2)') 10 format (1x, 500f12.4) 11 format (1x, 500f12.6) 15 format (1x, 3(f12.4,2x), f12.6) 16 format (1x, 2(f13.4,2x),2x, f13.4) 21 format (1x, 'top ', 1x, 4(f9.2,1x)) 22 format (1x, 'bottom', 1x, 4(f9.2,1x)) 23 format (1x, 'right ', 1x, 4(f9.2,1x)) 24 format (1x, 'left ', 1x, 4(f9.2,1x)) endif !if (btest(am0efl,1)) then !Erosion summary - total, salt/creep, susp, pm10 !(deposition values are positive - erosion values are negative) if (btest(am0efl,0)) then call caldatw (da, mo, yr) !get day, month and year write (UNIT=o_E_unit,FMT="(' ',i5,i3,i3,' ')",ADVANCE="NO") & & yr, mo, da write (UNIT=o_E_unit,FMT="(4(f12.6),' ')",ADVANCE="YES") & & asoilLossTot, (asoilLossTot-asoilLossSusp), asoilLossSusp, asoilLossPM10 endif ! end of plot section return end