Conference rusure::math

Let's try to fit a step voltage waveform to a normal distribution.  The 
purpose of this is to perform a chi-square test to determine the goodness of 
fit of the pulse to a Gaussian wave.  If it fits to .01 or .05, we can use 
shortcuts in correcting for bandwidth distortion.  If it doesn't fit, we may 
have to employ deconvolution techniques to recover a wave distorted by 
low-pass filtering of cables and connectors.

The waveform of the pulse input to the cable network looks like this:
(except that it's rounded and slightly exponential).

                /----------------
               /
______________/

The normal distribution is:
f(x)=1/(sigma * SQRT(2*pi)) * e**((-(x-u)**2)/(2*sigma**2))
where:
u = the mean and
sigma**2 = the variance.

The IEEE defines a Gaussian characteristic of a pulse such that it follows:
V = a*e**(-b(t-c)**2) 	(t>0) *******!

(We take V for voltage; the IEEE did not specify the parameter.)  To convert 
the standard equation to the IEEE we equate 
b=1/(2*sigma**2)   		Also,
a= 1/(sigma * SQRT(2*pi))   

sigma is in time units, so 'a' is in frequency squared units.

'c' is the time that the distribution would center around (the "average time")
if our pulse were a genuine normal distribution.  How should this point be 
selected?  At first I thought it would be the point at which the first 
derivative is zero in that part of the wave, that is, where it hits the
topline and is flat.  However, that does not allow a Gaussian fit.  The mean
time will intersect the pulse further down before it reaches the topline.
We also discard the pre-pulse epoch up to when the derivative is >0 and the
voltage rises above the baseline found by histogram, eliminating
pre-shoot (a slight drop just before rising).  The baseline can be
chosen by a histogram of the voltage.  

b also dependent on sigma, is equal to 1/(2*sigma**2).  Sigma in my case
was about 1/4 of the 20-80% risetime.

'a', the maximum voltage would be at the mean time of this bell.

The equation looks like this:
V = V(max) * e**(-((1/(2*(.25*tT)**2)*(t-(mean time))**2) 	(t>0)
tT is risetime, using 20%-80% proximal and distal points.

For a pulse from a fast logic chip through 6 feet of cable (I ran out of
time to attach the sampling head to the test boards, which was preferred,
    but the system seems to have a 54ps risetime),
I found the following equation, but visual inspection, because I couldn't
decide how to fit the data algorithmically:

V = .5 * EXP(-(1/(2*200E-12**2))*(T-.6E-9)**2)
and in my opinion it is a lousy fit, but I haven't done chi squared on it.
This equation forms a Gaussian curve under the pulse waveform that fits
the initial waveform rise up to about the 67% point and then suddenly deviates
to its peak, far short of the waveform topline, and returns to zero volts.

The shortcut correction for bandwidth distortion mentioned at the beginning
is used to find the input pulse's risetime duration from the output
risetime duration.  Given lossy cables and a sampling head of limited
bandwidth, the output risetime duration (tTo)will be longer, greater, than
the input risetime duration (tTi).
If the input pulse is Gaussian, then 
	tTo^2 = tTi^2  +  tTs^2
(I have no proof for this, do you?)
where tTs is the risetime of the cables and oscilloscope sampling head.
If the input isn't Gausssian in shape, this relationship doesn't hold
and the input risetime duration must be found by other methods.

If the Gaussian fit deviates from the waveform at the 67% level, then it
doesn't fit from 20% to 80%, the reference levels for transition time, 
which is the major application for the "Gaussian correction."

The result, although performed only graphically, seems to be that 
the waveform is not Gaussian and cannot use a Gaussian correction technique
to adjust for distortion due to low-pass filter functions of cables and
sampling heads.  Other techniques should be adapted, such as deconvolution.

\

P1q-?[a@@PPq??wSQPQSw??~!5?~??e!5Hq??e!5Hq???@@~@@???wSQPQSw??~ACGO_~??[a@@PPq??~HHHXhE??wSQPQSw??~!5HE??~!5G~!255?!216?CA~!5?[a@@@a[??!7G??wSQPQSw??~!5?~??[a@@PPq??!7G??u!5Hu??e!4He[!12?CA~!5?aPPHHIC!4?ee!5?nH!4Dx??OWSQ~OO!11?~!5HE??_!5Sw??W!4cwE??w!5SW
!12?CA~-???!4@???@!5?@???!5@!4?!5@!4?!5@!4?!5@???@!5?@??@!5?@!4?!4@???@!5?@??@!5?@??@!8?@!5?@!255?!216?!5@!5?@@@!13?@!5?@???!5@!5?!4@!13?!5@!4?!4@!14?!5@???!7@!4?@@!6?!5@!7?@!13?@!9?!6@??A!4CB!4?!6@!12?!5@--!185?oG!4CG???!5O_??o!4?o???_!5O???_!5O!5?Os!6
?!5O_??o_!4O_!11?{!4cCC!4?Os!5?OO{OO!13?oO_?_O_???!5O_??O_??_O???!7O!11?{c!4Sc!11?oO_?_O_??_!5O_???!5O_??o_!4O_??!7O!11?oG!4cG??o_!4O_??_!5O!12?_!5O!5?Os!5?_!4O_W??oO_?_O_???!5O_??!7O!11?GCCccgO??o_!4O_??_!5O-!185?@ACCDDB??A!5DF??B!4CBC??C!5DA??C!5DA???
CCFCC???A!5DF??F!5?F!11?F!9?CCFCC!5?BCCA!12?F??B??F??A!5DF??CA@@AC???!7@!4?EE!5?A!5CB!11?F??B??F??B!6D??A!5DF??F!5?F??!7@!4?EE!5?@B!4CB??F!5?F??C!5DA!11?C!5DA???CCFCC???H!4QN???F??B??F??A!5DF??!7@!4?EE!5?EDD!4C??F!5?F??C!5DA---!20?_OOO_!12?_!5O_-!19?FGO
OOGF!4?WW!5?L!5QL?!6A}!6A-!51?~-!51?~!255?!91?!19_o!19OW!17GK!35C!79ABB!7A@BBB!4AB!63@-!51?~!255?!15?!8_!9OW!12G!11CAEE!15AB!19@-!51?~!232?!4_o!5OW!5GK!5CE!7A!9@-!51?~!206?___o!4OWGGGK!4CAAAB!6@-!19?[a@@@a[!4?__!5?[q!4Hq?!6G~!6G!179?___oOOOWGGCCCEAAAB@@
@-!21?@@@!6?@@!7?!4@!8?~!169?__oOOWGKCCCEAB@@@-?_!5?_!43?~!154?___OOWGGCCEAB@@@-?@EW_WE@!43?~!137?___ooowwwWWW[[YYzx!5o__-??!5_!44?~!127?_oow{[mmRZJNDDEB@!14?@@BBBFM][wgoo_-?F!5GF!43?~!120?_oww[}vZHDCEB@!35?@@FJMK[wO_!255?!221?!19_???__???!19_!12?O__OOO
_??_!5O_???!5O_!4?C{-???A}!14?oGCCCGo!12?_OG{???!6O~!6O!108?_o{y]^LNCB@@!50?@BFFM]wo_!255?!255?!19?F!7?B!6D??A!5DF!4?CFC-???aBA!13?@ACCCA@!4?EE!5?!4@F@@!7?~!109?o{yMNFB@!65?@AAFN{go_!255?!209?!19W???GG???!19W!12?GG}HH@A!4?C|!5?CC~CC_-?AA^aaO!44?~!103?_o
{]NFB@!77?@BFNY{o_!255?!255?!9?@!7?!5@!6?@@-??!5_!44?~!96?_ow{k]NFA@!87?@AAF^sw_-?H!5IC!43?~!90?_oW{]yNFDB!99?BDDM}go-!51?~!83?__OWi}uTMNB@!109?BDFMY{o_-!20?_OOO_!12?_!4O_!8?~!77?__OWKtZ{]NDB!121?BENM{go_-!19?FGOOOGF!4?WW!5?WSSQQQP?!6@~!6@!65?_oOGIfxk[M
NB@!131?@BFF]sw_-!51?~!64?_oOOKerx[MNFB@!143?BDFEM[wo_-!51?~!58?OwgksQy|MFFBB!156?@AAFN]kwoo_-!51?~!49?_oog{}z\TNMFB@@!172?@FFMI[{woO__-!51?~!38?__ooww{[]MNNFBA@@!189?@@BBFNMK[[wwooo__-!21?___!15?___!9?~!11?!8_!5o!5WKKKMMEEFFBB@@@!81?_!124?_!7?@@BBBFEEM
MKK[[WWw!7o!6_!87?_!124?_!124?_-!19?MP___PM!4?oo!5?MP___PM?!6C~EE!14F!5D!103C~!124C~!29C!6D!13F!25E!51C~!124C~!124C~-!47?wCAAACw!109?]Q!4Iq??wCAAACw??wCAAACw!97?GC}!5?wCAAACw??wCAAACw??wCAAACw!92?GC}!5?]Q!4Iq??wCAAACw??wCAAACw!91?CaaQQSG??wCAAACw??wCAAA
Cw??wCAAACw!91?CaaQQSG??]Q!4Iq??wCAAACw??wCAAACw-!48?@AAA@!110?@!5A@???@AAA@!4?@AAA@!98?AABAA!4?@AAA@!4?@AAA@!4?@AAA@!93?AABAA???@!5A@???@AAA@!4?@AAA@!92?B!6A???@AAA@!4?@AAA@!4?@AAA@!92?B!6A??@!5A@???@AAA@!4?@AAA@--!255?!69?GGGwGGG!4?_g!5?__???_!4?!5_!1
4?_OG!4?_?!4_!4?!5_!5?GO_-!255?!72?N!6?GGNGG???N?@E@?N??F!5IJ!13?BCG!4?~B!4CB??H!5IC!4?GCB-\

    Prove the following theorem in statistics:
    
    Given the normal distribution 

    s=sigma
        f(x)= 1/(s*sqrt(2*pi)) * exp(-(x-u)^2 / (2 * s^2))
    
    Given the convolution integral:
    
    b(t) = [integral from 0 to t] p(t-tau)*g(tau) d-tau
    
    Prove that the convolution of the normal distribution with itself
    results in a normal distribution with a standard deviation s
    that is SQRT(2)*s of the original.
    
    I spent four days trying to look it up and solve it, but I couldn't
    even integrate the normal distribution.  Every book says,
    use the integral of the normal distribution to find the percentage
    below a certain x, but they never give you the integral.

    I don't think I solved the convolution integral either.
    I am convinced this theorem is right because I did it numerically
    and graphically.  hey, I'm a H/W engineer, OK? 8-)
    This is important for my daily work because there is misunderstanding
    about how filter function convolutions effect imperfect step waveforms
    and measurements of high-speed logic integrated circuits, but don't
    let that intrude on your fun.
    Tom

! GaussianConvolution.bas
! Tom Janzen 8-15-89
! convolutes a normal distribution with itself and display in ReGIS
MaxScreenX=768 
DIM WaveForm(MaxScreenX),Funct(MaxScreenX),convolved(MaxScreenX)
MaxScreenY=480

FOR WavePoint%=1 TO MaxScreenX
  s=(1/sqrt(2))*MaxScreenX/6
  u=MaxScreenX/2
  x=WavePoint%
  !f=1/(s*sqrt(2*pi)) * exp(-(x-u)**2/(2*(s**2)))
  f=exp(-(x-u)**2/(2*(s**2)))
  WaveForm(WavePoint%)=f
  Funct(WavePoint%)=f
  convolved(WavePoint%)=0
NEXT WavePoint%

print chr$(27);"[2J"
print chr$(27);"P1p"
print "S(E)";

startover:
GOSUB ShowWave  
GOSUB ShowFunction
GOSUB ConvolveWave
GOSUB ShowConvolved  
goto cleanup

ShowWave:
  FOR n=0 TO MaxScreenX-1
  screenx$ = str$(integer(int(n)))
  screeny$ = str$((MaxScreenY*.3) - integer(int((Waveform(n)) * MaxScreenY/4)))
  if n = 0 then print "p";
  END IF
  if n=1 then print "v";
  END IF
  print "[";screenx$;",";screeny$;"]";
  oldx$=screenx$
  oldy$=screeny$
next n
RETURN

ShowFunction:
  FOR n=0 TO MaxScreenX-1
  screenx$ = str$(integer(int(n)))
  screeny$ = str$((MaxScreenY*.6) - integer(int((Funct(n)) * MaxScreenY/4)))
  if n = 0 then print "p";
  END IF
  if n=1 then print "v";
  END IF
  print "[";screenx$;",";screeny$;"]";
  oldx$=screenx$
  oldy$=screeny$
next n

RETURN

ShowConvolved:
  FOR n=0 TO MaxScreenX-1
  screenx$ = str$(integer(int(n)))
  screeny$ = str$((MaxScreenY*.9) - integer(int((Convolved(n)) * MaxScreenY/4)))
  if n = 0 then print "p";
  END IF
  if n=1 then print "v";
  END IF
  print "[";screenx$;",";screeny$;"]";
  oldx$=screenx$
  oldy$=screeny$
next n

RETURN

ConvolveWave:
MaxCon=0
FOR WavePointer=0 TO MaxScreenX STEP 1
  TempWavePoint=0
  FOR FunctionPointer=0 TO WavePointer STEP 1
    TempWavePoint=(Funct(FunctionPointer))*(WaveForm(WavePointer-FunctionPointer))+TempWavePoint
  NEXT FunctionPointer
  convolved(WavePointer)=TempWavePoint
  IF TempWavePoint>MaxCon THEN MaxCon=TempWavePoint
  END IF
NEXT WavePointer
Expand=1/MaxCon 
FOR WavePointer=1 TO MaxScreenX STEP 1
  convolved(WavePointer)=convolved(WavePointer)*Expand
NEXT WavePointer
RETURN

cleanup:

print "p[10,10]T'Waveform'";
print "p[10,150]T'Function'";
print "p[10,300]T'Convolution of Waveform and Function'";
PRINT chr$(27);"\"

end

       <<< Note 1102.3 by ANT::JANZEN "cf. ANT::CIRCUITS,ANT::UWAVES" >>>
                        -< prove a Statistics theorem >-


>>     Prove the following theorem in statistics:
    
>>     Given the normal distribution 

>>     s=sigma
>>         f(x)= 1/(s*sqrt(2*pi)) * exp(-(x-u)^2 / (2 * s^2))
    
>>     Given the convolution integral:
    
>>     b(t) = [integral from 0 to t] p(t-tau)*g(tau) d-tau
    
>>     Prove that the convolution of the normal distribution with itself
>>     results in a normal distribution with a standard deviation s
>>     that is SQRT(2)*s of the original.
    
    The truth of the theorem follows from two statements that I know to be
    true (I can't look up the proof of the individual theorems that I quote,
    because all my math/statistics books are in a moving van somewhere between
    Atlanta and Colorado Springs):
    
    1.  The sum of two random variables with normal distributions is a
        random variable with a normal distribution.  If the means of the
        random variables are m1 and m2 and the standard deviations are s1
        and s2, then the mean of the sum is m1+m1 and the standard
        deviation of the sum is sqrt(s1*s1+s2*s2).
    
    2.  If two random variables have the frequency functions f1 and f2,
        then the sum of these random variables has a frequency function
        and that frequency function is the convolution of f1 and f2.
    
>>     I spent four days trying to look it up and solve it, but I couldn't
>>     even integrate the normal distribution.  Every book says,
>>     use the integral of the normal distribution to find the percentage
>>     below a certain x, but they never give you the integral.
    
    There is a very good reason why the books never quote the indefinite
    integral of the normal distribution function: it can not be given as a
    finite sum of elementary functions.  In practice, one evaluates the
    integral numerically.

>>     I don't think I solved the convolution integral either.
    
    I haven't looked closely at it, but it looks likely that the
    convolution integral will wind up very similar in form to the integral
    of the normal-distribution frequency function; my advice would be to
    not even try to find an indefinite integral.
    
>>     I am convinced this theorem is right because I did it numerically
>>     and graphically.  hey, I'm a H/W engineer, OK? 8-)
    
    Nothing wrong with being a H/W engineer!
    
>>     This is important for my daily work because there is misunderstanding
>>     about how filter function convolutions effect imperfect step waveforms
>>     and measurements of high-speed logic integrated circuits, but don't
>>     let that intrude on your fun.
>>     Tom

    If you need more help, drop me a line or give me a call.

        I sketched the proof out on paper, and will try to
        enter it sometime this weekend.  It takes time to
        enter anything involving nontrivial notation.
        
        The only step not completely justified in the "sketch"
        is a change in the order of taking the derivative of an
        integral:  d/dy int f(y,t)dt <--> int f (y,t) dt
                                               y
        Dan

	A real-valued random variable X can be specified by its
	cumulative distribution function Prob(X <= x) = F(x).
	It's density function is then the derivative of that,
	f(x) = F'(x), so that Prob(X <= x) = F(x) = the integral
	from -oo to x of f(t) dt.
							      2
	For the normal distribution with mean u and variance s
	you have f(x) = (1/(s sqrt(2 pi))) exp(- (x - u)^2 / (2s^2))
	and F is the integral of that as before, it has no closed
	form solution.  I think the usual designation for the
	cumulative distribution function of a normal distribution is
	the "erf" function, where F(x) is either erf( (x-u)/s ) or
	that plus 1/2.

	If X and Y are random variables with cdf's F and G and
	density functions f and g, then the cdf for the derived
	random variable Z = X + Y is given by a double integral
	over the region of the x,y plane that is left/below the
	line x + y = z ...

	Prob(Z <= z) = H(z) =
	int(x from -oo to +oo) f(x) dx int(y from -oo to z-x) g(y) dy
	= int(x from -oo to +oo) f(x) G(z-x) dx

	Now you get the density function h(z) for Z by taking the
	derivative of H(z) with respect to z.  If you move the
	differentiation inside the integral (hey, we're all
	friends here, no one should mind! :-) ) then you get

	h(z) = (d/dz) int(x from -oo to +oo) f(x) G(z-x) dx
	     = int(x from -oo to +oo) f(x) g(z-x) dx

	and this is the convolution of the density functions f(x)
	and g(y) of the random variables X and Y.  So if X and Y
	are normally distributed with mean and standard deviation
	u and s for X and v and t for Y, then

		h(z) = int(x from -oo to +oo) f(x) g(z - x) dx
		f(x) = (1/(s sqrt(2 pi))) exp( - (x - u) / (2s^2) )
		g(y) = (1/(t sqrt(2 pi))) exp( - (y - v) / (2t^2) )

	Now this is a mess to try to write down on a VT220, but
	here goes.  Any terms that do not have an x in them can
	be brought outside of the integral.  The product of the
	two exponentials becomes a single exponential of the sum
	of the exponents.  This becomes two exponentials again,
	one inside the integral with exponent the sum of the
	terms involving x, one outside the integral with exponent
	the terms not involving x.

	Now, the constant factors do not involve x and result in
	a factor of 1/(2 pi s t) outside of the integral.  The
	combined exponent of the exp(...) terms is minus

	       2	      2
	(x - u)    (z - x - v)
	-------- + ------------ = minus
	  2s^2         2t^2

	t^2 x^2 - 2 t^2 x u + t^2 u^2 + s^2 (z-v)^2 - 2 s^2 (z-v) x + s^2 x^2
	---------------------------------------------------------------------
				2 s^2 t^2

	The part that can be moved outside the integral is exp of
	minus (t^2 u^2 + s^2 (z-v)^2) / (2 s^2 t^2).  The part that
	stays inside the integral is exp of minus

	((s^2 + t^2) x^2 - 2 (s^2 (z-v) + t^2 u) x) / (2 s^2 t^2)

	Divide through numerator and denominator by (s^2 + t^2) and
	let a = (s^2 (z-v) + t^2 u)/(s^2 + t^2),
	let b^2 = (s^2 t^2)/(s^2 + t^2), to get

		exp( - (x^2 - 2ax) / (2b^2))

	inside the integral.  This can be solved as follows:
	complete the square in the numerator by multiplying
	inside the integral by exp( -a^2 / (2b^2)) and outside
	the integral by exp(a^2 / (2b^2)).  Then the integral is
	int(x from -oo to +oo) exp (- (x-a)^2 / (2b^2)) which
	is b sqrt(2 pi) by comparing it with the form of the
	normal density function f(x) above and noting that
	int(x from -oo to +oo) f(x) dx = 1 because all cdf's
	approach 0 towards -oo and approach 1 towards +oo.

	Combining all of the terms, we have the result that h(z)
	is the product of 1/(2 pi s t), exp of minus
	(t^2 u^2 + s^2 (z-v)^2) / (2 s^2 t^2), b sqrt(2 pi), and
	exp(a^2 / (2b^2)).  This should work out to a normal with
	mean u+v and variance s^2 + t^2.  I don't feel like
	showing that, so I will let s=t (all that was requested in .3)
	and simplify.  Then a = (z-v+u)/2, b^2 = (s^2)/2, st = s^2,
	(t^2 u^2 + s^2 (z-v)^2) / (2 s^2 t^2) =
	(u^2 + z^2 - 2zv + v^2) / (2 s^2), and the result is

	(1/(2 pi s^2)) (s/sqrt(2)) sqrt(2 pi) times exp of
	-(u^2 + z^2 - 2zv + v^2 - 2 a^2) / (2 s^2)

	= (1 / (s sqrt(2) sqrt(2 pi))) exp( -(u^2 + z^2 - 2zv + v^2
	-(1/2) z^2 + zv - zu + uv -(1/2) - (1/2) u^2 - (1/2) z^2) /
	(2 s^2))

	= (1 / ( (s sqrt(2)) sqrt(2 pi) )) exp( - (z^2 - 2zv - 2zu
	+ u^2 + 2uv + v^2) / (2 (s sqrt(2))^2))

	Let S = s sqrt(2) and this is

		(1/(S sqrt(2))) exp( - (z - (u+v))^2 / (2S^2))

	which is exactly a normal distribution with mean u+v and
	variance S^2 = 2s^2, which was the result requested in .3.

	Dan

        re .0
        
>>    I can't find an algorithm for a gaussian distribution fit in my
>>    books.  How do I fit data and do a chi-square test for the equation
>>    
>>    V= A * e ** (-b(t-c)**2)
>>    ?
        
        Maybe I am missing something, but the usual way to fit a
        Gaussian to n data points is to take the mean u of the
        Gaussian to be m = (x1 + ... + xn) / n, and the variance s^2
        of the Gaussian to be ((x1 - m)^2 + ... (xn - m)^2) / n
        [or perhaps divided by (n-1) instead of by n].  Which of
        V, A, b, and c do you have data for, and which are you
        trying to estimate?  Of the ones with data which is/are
        "independent" and which "dependent"?
        
        Dan

    Hi
    Thanks for the convolution information.
    Anyone who wants the vms basic regis software can ask me. It grahps
    two gaussian curves and convolutes them and draws the result.
    
    What I am really trying to fit is a real step waveform, which is
    kind of curvy of course.  I want to fit the transition epoch to
    the left part of a gaussian distribution.
                   /----------------
     -------------/
                  ___
      -----------/   \_________________

    So I know the discrete numerical 
    density function, V(t) and want to fit it to
    V(t)=A * e**(-b*(t-c)**2))

    Thanks
    Tom

    This was also part of my reading on convolution and deconvolution.
    VAX BASIC programs for both operations are in ant::circuits 203
    Tom