c $Header: /weru/cvs/weps/weps.src/asd/m2asd.for,v 1.4 1999-04-22 19:02:30 wjr Exp $ c c NOTES: c c We will probably want to rewrite what this subroutine does into c several smaller routines for both speed reasons and potentially c modularity. c c The log(md) that we need could be moved to the c the initialization routine and access the log(md) rather than c accessing (md) and computing log(md) multiple times (SPEED SAVINGS). c c We also may want to may want to separate out the looping among c soil layers so that it can be done at a higher level (may make c code more modular - do only one thing extremely well concept). c This should be discussed as to whether this would be beneficial c in the long run. c c Tue Apr 6 14:15:48 CDT 1999 - LEW c ----------------------------------------------------------------- c This routine was simplified and recoded. c It now allows for the sieve cut sizes to lie outside the c range specified by "mnot" and "minf" by checking for this c situation and only using the sieve cuts between "mnot" and "minf" c (this only applies to the pertinent modified log-normal cases). c c Note that: c a) the sieve cut size array, "mdia" must consist of 2 or more sizes c and contain values which increase in size, c b) "mnot" must be greater than or equal to zero, c c) "mnot" must be greater than "minf" (with at least two sieve cut c sizes between them), c d) and the mass fractions, "mf" cannot be less than zero. c These conditions are NOT checked within this code. c c Note also that the return values "gmd" and "gsd" are the c geometric mean and geometric standard deviation of the c "transformed" parameters, based upon the specific "logcas" c used and NOT always the geometric mean and standard deviation c of the aggregate sizes. c ----------------------------------------------------------------- subroutine m2asd i (mf, nlay, i mnot, minf, o gmd, gsd) *$noereference include 'p1werm.inc' include 'manage/asd.inc' *$reference c + + + PURPOSE + + + c This subroutine performs the inverse of subroutine asd2m. c m2asd computes the geometric mean & standard deviation for the c lognormal representation of the soil aggregate size distribution c from mf(i,j). c The routine decides which lognormal case to apply based on the c value of logcas: c c logcas = 0 --> "normal" lognormal case (mnot = 0, minf = infinity) c logcas = 1 --> "abnormal" lognormal case (mnot != 0, minf = infinity) c logcas = 2 --> "abnormal" lognormal case (mnot = 0, minf != infinity) c logcas = 3 --> "abnormal" lognormal case (mnot != 0, minf != infinity) c c + + + KEYWORDS + + + c aggregate size distribution, asd, sieves, mass fractions c c + + + ARGUMENT DECLARATIONS + + + real mf(msieve+1, mnsz) integer nlay real mnot(mnsz), minf(mnsz) real gmd(mnsz), gsd(mnsz) c c c + + + ARGUMENT DEFINITIONS + + + c mf - mass fractions of aggregates within sieve cuts c (sum of all mass fractions are expected to = 1.0) c nlay - number of soil layers used c mnot - minimum size aggregate (assumed value is known) c minf - maximum size aggregate (assumed value is known) c gmd - geometric mean diameter of aggregate size distribution c (or transformed asd for "modified" lognormal cases) c gsd - geometric standard deviation of aggregate size distribution c (or transformed asd for "modified" lognormal cases) c c + + + ACCESSED COMMON BLOCK VARIABLE DEFINITIONS + + + c c nsieve - number of sieves used c mdia - geometric mean dia. for each sieve cut c mnsize - minimum (imaginary) sieve size to use for computing c lower sieve cut geometric mean diameter c mxsize - maximum (imaginary) sieve size to use for computing c upper sieve cut geometric mean diameter c logcas - flag to represent which lognormal case to apply c c + + + PARAMETERS + + + c c + + + LOCAL VARIABLES + + + c real tmd(msieve+1) real alpha, beta integer i, j integer istart, istop c c + + + LOCAL VARIABLE DEFINITIONS + + + c c tmd - transformed md (later log(tmd)) c alpha - internal summation variable c beta - internal summation variable c i - loop variable for sieve diameters c istart - loop start variable for sieve diameters c istop - loop stop variable for sieve diameters c j - loop variable for soil layers c c + + + END SPECIFICATIONS + + + c c for each soil layer do 200 j=1,nlay c initialize accumulators alpha = 0.0 beta = 0.0 istart = 1 istop = nsieve + 1 c check if sieve cut fractions are between mnot and minf do 10 i=nsieve+1, 1 if (mdia(i) .ge. mnot(j)) then istart = i end if 10 continue if (logcas .ge. 2) then do 20 i=1, nsieve+1 if (mdia(i) .le. minf(j)) then istop = i end if 20 continue else istop = nsieve + 1 end if c do transformations for "modified" log-normal cases do 30 i= istart, istop if (logcas .ge. 2) then tmd(i) = (mdia(i)-mnot(j))*(minf(j)-mnot(j))/ & (minf(j)-mdia(i)) else tmd(i) = mdia(i)-mnot(j) end if c now compute the log of the gmd dia tmd(i) = log(tmd(i)) c sum diameters & their squares, over all aggregate sizes alpha = alpha + (mf(i,j)*tmd(i)) beta = beta + (mf(i,j)*tmd(i)*tmd(i)) 30 continue c compute geometric mean and standard deviation gmd(j) = exp(alpha) C *** debugging fix if (beta-alpha**2.lt.0) then gsd(j) = 1 write(*,*) 'beta, alpha ', beta, alpha C *** eodf else gsd(j) = exp(sqrt(beta-alpha*alpha)) endif 200 continue return end c c$Log: not supported by cvs2svn $ cRevision 1.3 1999/04/06 19:40:29 wagner cSimplified and recoded this routine. cIt now checks for and allows sieve fractions to be outside cthe range from "mnot" to "minf" for the "modified" log-normal ccases. c cRevision 1.1.1.1 1999/03/12 17:05:17 wagner cBaseline version of WEPS with Bill Rust's modifications c c Revision 1.2 1995/09/13 15:49:44 wagner c Necessary changes made to allow FORTRAN src files (*.for) to use the c extended FORTRAN include statement rather than the MICROSOFT $INCLUDE c directive as previously used. This is required to allow use of other c FORTRAN compilers. c c Changes have been made to the prologue.mk, epilogue.mk, and the Unix c master startup.mk files as well as the src files. c c Revision 1.1.1.1 1995/01/18 04:19:56 wagner c Initial checkin c c Revision 1.4 1994/09/01 22:18:54 jt c checking for floating point errors? - LEW c c Revision 1.3 1992/10/10 21:44:14 wagner c Changed names appropriate for submodel name change c from TILLAGE to MANAGEMENT. c c Revision 1.2 1992/04/16 21:41:37 wagner c Uses common memory now. c c Revision 1.1 1992/04/16 13:29:01 wagner c Initial revision c