[python] Re: Bicycling physics

Hi Ray,

when I ran your program I got a little different output than the one you copied pasted in your mail.

I had also to change the for loop in the gtrail (although it prints out trails going from 10 -> -25cm, the program uses -40 -> -48cm). With 10->-25 the values printed out get quite wierd (no self centering effect):
steerAngle 20.0
pivot                           ground trail                      
        10      5       0       -5      -10     -15     -20     -25

40      -4.077  -3.532  -2.98   -2.421  -1.855  -1.281  -0.699  -0.111 
41      -4.109  -3.562  -3.007  -2.444  -1.873  -1.295  -0.708  -0.114 
42      -4.139  -3.589  -3.031  -2.466  -1.892  -1.31   -0.719  -0.12 
43      -4.165  -3.614  -3.055  -2.487  -1.91   -1.325  -0.731  -0.129 
44      -4.189  -3.637  -3.077  -2.507  -1.929  -1.342  -0.745  -0.14 
45      -4.209  -3.658  -3.097  -2.527  -1.947  -1.359  -0.761  -0.153 
46      -4.226  -3.675  -3.115  -2.545  -1.966  -1.376  -0.777  -0.168 
47      -4.24   -3.69   -3.131  -2.562  -1.983  -1.394  -0.795  -0.185 
48      -4.249  -3.702  -3.145  -2.578  -2.0    -1.412  -0.813  -0.204 
49      -4.255  -3.711  -3.156  -2.591  -2.016  -1.43   -0.833  -0.224 
50      -4.257  -3.716  -3.165  -2.604  -2.031  -1.447  -0.852  -0.246 
51      -4.254  -3.718  -3.171  -2.614  -2.045  -1.465  -0.873  -0.269 
52      -4.248  -3.717  -3.175  -2.622  -2.057  -1.481  -0.893  -0.293 
53      -4.236  -3.711  -3.175  -2.627  -2.068  -1.497  -0.914  -0.319 
54      -4.221  -3.702  -3.172  -2.631  -2.077  -1.512  -0.935  -0.344 
55      -4.2    -3.689  -3.166  -2.631  -2.085  -1.526  -0.955  -0.371 
56      -4.175  -3.671  -3.156  -2.629  -2.09   -1.539  -0.975  -0.398 
57      -4.145  -3.649  -3.142  -2.624  -2.093  -1.55   -0.994  -0.424 
58      -4.109  -3.623  -3.125  -2.616  -2.094  -1.559  -1.012  -0.451 
59      -4.069  -3.592  -3.104  -2.604  -2.092  -1.567  -1.029  -0.478 
60      -4.023  -3.557  -3.079  -2.589  -2.087  -1.572  -1.045  -0.504 
61      -3.972  -3.517  -3.05   -2.571  -2.08   -1.576  -1.059  -0.529 
62      -3.916  -3.471  -3.016  -2.548  -2.069  -1.577  -1.072  -0.554 
63      -3.854  -3.421  -2.978  -2.522  -2.055  -1.575  -1.083  -0.577 
64      -3.786  -3.366  -2.935  -2.493  -2.038  -1.571  -1.092  -0.599 
65      -3.713  -3.306  -2.888  -2.458  -2.017  -1.564  -1.098  -0.619 
66      -3.634  -3.241  -2.836  -2.42   -1.993  -1.554  -1.102  -0.638 
67      -3.55   -3.17   -2.779  -2.378  -1.965  -1.54   -1.103  -0.654 
68      -3.46   -3.094  -2.717  -2.33   -1.932  -1.523  -1.102  -0.668 
69      -3.364  -3.012  -2.651  -2.279  -1.896  -1.503  -1.097  -0.68 
70      -3.262  -2.925  -2.579  -2.222  -1.856  -1.478  -1.089  -0.689 
71      -3.154  -2.832  -2.502  -2.161  -1.811  -1.45   -1.078  -0.695 
72      -3.04   -2.734  -2.419  -2.095  -1.762  -1.418  -1.063  -0.698 
73      -2.921  -2.631  -2.332  -2.024  -1.708  -1.381  -1.045  -0.698 
74      -2.795  -2.521  -2.239  -1.948  -1.649  -1.34   -1.022  -0.694 
75      -2.664  -2.406  -2.141  -1.867  -1.586  -1.295  -0.995  -0.686 
76      -2.527  -2.286  -2.037  -1.781  -1.517  -1.245  -0.964  -0.675 
77      -2.384  -2.16   -1.928  -1.69   -1.444  -1.191  -0.929  -0.659 
78      -2.235  -2.028  -1.814  -1.593  -1.366  -1.131  -0.889  -0.639 
79      -2.081  -1.89   -1.694  -1.491  -1.282  -1.067  -0.845  -0.615
To be able to compare, I adapted my program (octave). Here is the output of my program:
octave> flevo2
aa	10.0	5.0	0.0	-5.0	-10.0	-15.0	-20.0	-25.0
40.0	-0.26	-0.12	0.01	0.15	0.29	0.42	0.56	0.69
41.0	-0.26	-0.12	0.01	0.15	0.29	0.42	0.56	0.70
42.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.56	0.70
43.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.57	0.70
44.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.57	0.71
45.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.57	0.71
46.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.57	0.71
47.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.57	0.71
48.0	-0.26	-0.12	0.01	0.15	0.29	0.43	0.56	0.70
49.0	-0.26	-0.12	0.01	0.15	0.29	0.42	0.56	0.70
50.0	-0.26	-0.12	0.01	0.15	0.29	0.42	0.56	0.70
51.0	-0.25	-0.12	0.01	0.15	0.28	0.42	0.56	0.69
52.0	-0.25	-0.12	0.01	0.15	0.28	0.42	0.55	0.69
53.0	-0.25	-0.12	0.01	0.15	0.28	0.41	0.55	0.68
54.0	-0.25	-0.12	0.01	0.14	0.28	0.41	0.54	0.67
55.0	-0.25	-0.12	0.01	0.14	0.27	0.40	0.53	0.67
56.0	-0.24	-0.12	0.01	0.14	0.27	0.40	0.53	0.66
57.0	-0.24	-0.11	0.01	0.14	0.26	0.39	0.52	0.65
58.0	-0.24	-0.11	0.01	0.14	0.26	0.39	0.51	0.64
59.0	-0.23	-0.11	0.01	0.13	0.26	0.38	0.50	0.63
60.0	-0.23	-0.11	0.01	0.13	0.25	0.37	0.49	0.62
61.0	-0.22	-0.11	0.01	0.13	0.25	0.36	0.48	0.60
62.0	-0.22	-0.10	0.01	0.12	0.24	0.36	0.47	0.59
63.0	-0.21	-0.10	0.01	0.12	0.23	0.35	0.46	0.58
64.0	-0.21	-0.10	0.01	0.12	0.23	0.34	0.45	0.56
65.0	-0.20	-0.10	0.01	0.11	0.22	0.33	0.44	0.54
66.0	-0.20	-0.09	0.01	0.11	0.21	0.32	0.42	0.53
67.0	-0.19	-0.09	0.01	0.11	0.21	0.31	0.41	0.51
68.0	-0.18	-0.09	0.01	0.10	0.20	0.30	0.40	0.49
69.0	-0.18	-0.09	0.01	0.10	0.19	0.29	0.38	0.48
70.0	-0.17	-0.08	0.01	0.09	0.18	0.27	0.37	0.46
71.0	-0.16	-0.08	0.01	0.09	0.18	0.26	0.35	0.44
72.0	-0.16	-0.08	0.00	0.09	0.17	0.25	0.33	0.42
73.0	-0.15	-0.07	0.00	0.08	0.16	0.24	0.32	0.40
74.0	-0.14	-0.07	0.00	0.08	0.15	0.23	0.30	0.38
75.0	-0.13	-0.07	0.00	0.07	0.14	0.21	0.28	0.36
76.0	-0.13	-0.06	0.00	0.07	0.13	0.20	0.27	0.33
77.0	-0.12	-0.06	0.00	0.06	0.13	0.19	0.25	0.31
78.0	-0.11	-0.05	0.00	0.06	0.12	0.17	0.23	0.29
79.0	-0.10	-0.05	0.00	0.05	0.11	0.16	0.21	0.27

I think the calculation of the frontAxle* values is too simplified. msssHT.py uses a linear formula for that, but that makes a rather large error. Instead of:

    k=steerAngle*2/pi  

I would propose

    k=1 - cos(steerAngle)

This has considarable effect. For example at 20degrees k is only 0.0603 instead of 0.222


The second problem is in the steerDropAngle. The correct formula is rather complicated. I switch to a matrix calculation here.

if you define the unrotated front wheel as [xt,yt,zt] =

    [R*cos(t), 0, R*sin(t)]   where t in [0, 2*pi[ and R is the front wheel radius

then the formula for the front wheel rotated in the new steering axle point would be (choke!):

    [cos(HA)*cos(HA)*R*cos(t) + cos(SA)*(1-cos(HA)*cos(HA))*R*cos(t) -sin(HA)*cos(HA)*R*sin(t) + cos(SA)*sin(HA)*cos(HA)*R*sin(t),
     -sin(SA)*sin(HA)*R*cos(t)- sin(SA)*cos(HA)*R*sin(t),
    -cos(HA)*sin(HA)*R*cos(t) + cos(SA)*cos(HA)*sin(HA)*R*cos(t) + sin(HA)*sin(HA)*R*sin(t) + cos(SA)*(1-sin(HA)*sin(HA))*R*sin(t)]

where HA= headangle, SA is steerAngle

This formula follows from substituting the formulas at:
    http://www.euclideanspace.com/maths/algebra/matrix/orthogonal/rotation/

The minimum z would be at (follows out of the derivation of the z component in parameter t

    tg(t) = ((sin(HA)*sin(HA) + cos(SA)*(1-sin(HA)*sin(HA))) / (cos(SA)*cos(HA)*sin(HA) - cos(HA)*sin(HA)))

This t you need to put back into the largish formula above to find the Z/X value of the lowest point of the wheel. This values can then be used to calculate the steerDropAngle.

As a test, I adapted your program (molested it actaully), but there are still a difference between the results. The new results of your program are here:

steerAngle 20.0
pivot                           ground trail                      
        10      5       0       -5      -10     -15     -20     -25

40      -0.274  -0.137  0.0     0.138   0.276   0.415   0.555   0.694 
41      -0.275  -0.138  0.0     0.139   0.278   0.418   0.558   0.698 
42      -0.276  -0.138  0.0     0.139   0.279   0.419   0.56    0.701 
43      -0.277  -0.139  0.0     0.14    0.28    0.421   0.562   0.703 
44      -0.277  -0.139  0.0     0.14    0.28    0.421   0.563   0.705 
45      -0.277  -0.139  0.0     0.14    0.28    0.421   0.563   0.705 
46      -0.277  -0.139  0.0     0.14    0.28    0.421   0.562   0.705 
47      -0.276  -0.138  0.0     0.14    0.279   0.42    0.561   0.703 
48      -0.275  -0.138  0.0     0.139   0.278   0.419   0.559   0.701 
49      -0.274  -0.137  0.0     0.138   0.277   0.416   0.557   0.698 
50      -0.273  -0.137  -0.0    0.137   0.275   0.414   0.553   0.693 
51      -0.271  -0.136  -0.0    0.136   0.273   0.411   0.549   0.689 
52      -0.269  -0.135  -0.001  0.135   0.271   0.407   0.545   0.683 
53      -0.266  -0.134  -0.001  0.133   0.268   0.403   0.539   0.676 
54      -0.264  -0.133  -0.001  0.131   0.264   0.398   0.533   0.669 
55      -0.261  -0.132  -0.001  0.129   0.261   0.393   0.526   0.66 
56      -0.258  -0.13   -0.002  0.127   0.257   0.388   0.519   0.651 
57      -0.254  -0.129  -0.002  0.125   0.253   0.381   0.511   0.641 
58      -0.25   -0.127  -0.003  0.122   0.248   0.375   0.502   0.63 
59      -0.246  -0.125  -0.003  0.12    0.243   0.368   0.493   0.619 
60      -0.242  -0.123  -0.003  0.117   0.238   0.36    0.483   0.606 
61      -0.237  -0.121  -0.004  0.114   0.233   0.352   0.472   0.593 
62      -0.233  -0.119  -0.004  0.111   0.227   0.344   0.461   0.58 
63      -0.227  -0.116  -0.005  0.108   0.221   0.335   0.45    0.565 
64      -0.222  -0.114  -0.005  0.104   0.215   0.326   0.437   0.55 
65      -0.216  -0.111  -0.006  0.101   0.208   0.316   0.425   0.534 
66      -0.211  -0.109  -0.006  0.097   0.201   0.306   0.411   0.518 
67      -0.205  -0.106  -0.007  0.093   0.194   0.295   0.397   0.5 
68      -0.198  -0.103  -0.007  0.089   0.187   0.285   0.383   0.483 
69      -0.192  -0.1    -0.008  0.085   0.179   0.273   0.368   0.464 
70      -0.185  -0.097  -0.008  0.081   0.171   0.262   0.353   0.445
For lower pivot angles the results are now more or less the same as in my program, but at larger pivot angles there is a diviation. But at least the values that are output are sane now (negative trail give a rise at steering angle).

I haven't checked the reason for the deviation (and like i said it is not impossible that I am wrong. Calculating all this I find quite tricky).

Attached is the adapted program.  Tip: you might have a look at quaternions. They simplify working in three dimensional space a great deal. AFAIK, there is a quaternion package in the python programming language.

cheers
Dirk






mtb@xxxxxxx wrote:
Hi Dirk,

Attached is the code I referred to. It uses Jürgen's newest P* geometry
Attached also is a sketch I used as a guide.
    wheelBase = 72.5
    steerAngle = 20.  * math.pi / 180.0
    rearMassZ = 28.
    rearMassX = 67.
    frontWheelRadius = 33.
    rearWheelRadius = 33.

The code tries to calculate the rise of the rear mass using a methodology I saw in one of the PDFs.
Some assumptions:
-a wheelie is a positive angle, the rear wheel is fixed and the frame rotates around it
-the rear frame does not lean
-turnAnglePnt is imaginary point the axle rotates about with steer,
is fixed with respect to rear frame, and is the center of a torus
swept by the wheel. It is on the pivot axis perpendicular to the front hub.
-all measurements are in the x-z plane

In the first part, measure the x-z change of the front hub as the rear frame is fixed and the pivot rotates. The front axle moves along the x-z plane with steer, sweeping through the torus. Then calc the angle between initial 0 steer, and current steer, relative to rear axle height; add any existing angle between axles. (If the front wheel was a point, this would definitely be the frame rise angle.)

Now, rotate the steer without mechanical trail holding frame fixed, and calc the low z point of the tire which change due to wheel flop. Calc angle due to steer -
max fWheel drop/flop is proportional to steer Angle==90

Subtract the angle the wheel drops due to flop, from the angle the frame rises due to the whole wheel sweeping about the center of the torus
Calc the angle from the rear mass to the rear axle.
Then:
newrearMassZ = rearWheelRadius+
  rmass2Raxle*sin(rmass2RaxleAngle+frontAxleAngle-steerDropAngle)

If you see something I missed, please let me know! Each part checks out for the geometric extremes, steer=0,90,180, but the "peak rise" result differs from your reference frame vector method.

By the way, I saw a very good German bike site:
http://www.kreuzotter.de/english/elenk.htm
Maybe some of you know Walter Zorn?

Ray

sample rise output below.


steerAngle 20.0
                              ground trail
pivot   10      5       0       -5      -10     -15     -20     -25
angle
40      5.753   5.862   5.971   6.081   6.191   6.301   6.411   6.521
41      5.787   5.898   6.008   6.119   6.23    6.341   6.452   6.564
42      5.817   5.928   6.04    6.151   6.263   6.375   6.488   6.6
43      5.841   5.953   6.066   6.178   6.291   6.404   6.518   6.631
44      5.86    5.973   6.086   6.2     6.313   6.427   6.541   6.656
45      5.874   5.988   6.102   6.216   6.33    6.445   6.559   6.674
46      5.883   5.997   6.112   6.226   6.341   6.456   6.571   6.687
47      5.887   6.001   6.116   6.231   6.346   6.461   6.577   6.693
48      5.885   6.0     6.114   6.23    6.345   6.461   6.577   6.693
49      5.878   5.993   6.108   6.223   6.338   6.454   6.57    6.687
50      5.865   5.98    6.095   6.21    6.326   6.442   6.558   6.674
51      5.848   5.962   6.077   6.192   6.307   6.423   6.539   6.655
52      5.824   5.939   6.053   6.168   6.283   6.398   6.514   6.63
53      5.796   5.909   6.023   6.138   6.252   6.367   6.483   6.598
54      5.762   5.875   5.988   6.102   6.216   6.33    6.445   6.56
55      5.723   5.835   5.947   6.06    6.174   6.287   6.401   6.516
56      5.678   5.789   5.901   6.013   6.125   6.238   6.351   6.465
57      5.628   5.738   5.849   5.96    6.071   6.183   6.295   6.408
58      5.573   5.682   5.791   5.901   6.011   6.121   6.232   6.344
59      5.513   5.62    5.728   5.836   5.945   6.054   6.164   6.274
60      5.447   5.553   5.659   5.766   5.873   5.981   6.089   6.198
61      5.376   5.48    5.585   5.69    5.796   5.902   6.008   6.115
62      5.3     5.403   5.506   5.609   5.713   5.817   5.922   6.027
63      5.22    5.32    5.421   5.522   5.624   5.726   5.829   5.932
64      5.134   5.232   5.331   5.43    5.53    5.63    5.731   5.832
65      5.044   5.14    5.236   5.333   5.43    5.528   5.626   5.725
66      4.949   5.042   5.136   5.231   5.325   5.421   5.517   5.613
67      4.849   4.94    5.031   5.123   5.215   5.308   5.401   5.495
68      4.745   4.833   4.922   5.011   5.1     5.19    5.281   5.372
69      4.637   4.722   4.808   4.894   4.981   5.068   5.155   5.243
70      4.524   4.606   4.689   4.772   4.856   4.94    5.025   5.11
71      4.408   4.487   4.566   4.646   4.727   4.808   4.889   4.971
72      4.287   4.363   4.439   4.516   4.593   4.671   4.749   4.827
73      4.163   4.236   4.308   4.382   4.455   4.529   4.604   4.679
74      4.035   4.104   4.174   4.243   4.313   4.384   4.455   4.526
75      3.905   3.97    4.035   4.101   4.168   4.234   4.302   4.369
76      3.77    3.832   3.894   3.956   4.018   4.081   4.145   4.208
77      3.633   3.691   3.749   3.807   3.866   3.925   3.984   4.044
78      3.494   3.547   3.601   3.655   3.71    3.765   3.82    3.875
79      3.351   3.401   3.45    3.5     3.551   3.602   3.653   3.704




At 08:56 AM 3/10/2005 +0100, you wrote:
Hi Ray,

OK, can you sent it to me (or better to the group). This weekend I should be able to look at it. Do note that it is not impossible that the error is on my side ;-)

It has been very busy on the list!!! I have only had the chance to look at what has been going on, without having the time to delve into the discussion. But I am looking forward to catch up with you guys!

Good stuff!
Dirk

Ray Schümacher wrote:

Hello Dirk,

I had coded a short routine in Python to calculate rear mass rise, but I still get somewhat different results than yours. Might you have time to look at it and see if anything jumps out at you? I have an image as well, sketching out the geometry. I'll post them to the group if that would be helpful.

Cheers,
Ray


 



#!/usr/bin/env python
#Boa:PyApp:main

import math

modules ={}

def main():

    wheelBase = 72.5
    steerAngle = 20.  * math.pi / 180.0
    rearMassZ = 28.
    rearMassX = 67.
    frontWheelRadius = 33.
    rearWheelRadius = 33.


    ## calc angle from raxle to rmass, from horizontal
    rmass2RaxleAngle = math.asin((rearMassZ-rearWheelRadius)/rearMassX)
    ## calc dist from raxle to rmass
    rmass2Raxle = ((rearWheelRadius-rearMassZ)**2+rearMassX**2)**.5
    ## angle front to rear axle
    ## a wheelie is a pos+ angle
    angleBtwnAxles = math.atan((frontWheelRadius-rearWheelRadius)/wheelBase)
    
    print 'steerAngle', steerAngle/math.pi*180.0
    print 'pivot                           ground trail                      '
    print '     ', '\t10', '\t5', '\t0', '\t-5', '\t-10', '\t-15', '\t-20', '\t-25'
    for angle in range(40, 80, 1):        
        headAngle = angle  * math.pi / 180.0 # deg from horiz
        print '\n',int(round(headAngle*180/math.pi)), '   ',
        ## iterate over the ground trail
        for gtrail in range(-42, -50, -1):            
            mechTrail = frontWheelRadius*math.cos(headAngle) - gtrail*math.sin(headAngle) 
                
            ## rear wheel does not move or roll, frame rotates around the axle
            ## turnAnglePnt is imaginary point the axle rotates about with steer, 
            ## is fixed with respect to rear frame, and is the center of a torus
            ## swept by the wheel
            turnAnglePntZ = frontWheelRadius+(math.cos(math.pi-headAngle)*mechTrail)#
            turnAnglePntX = wheelBase-(math.sin(math.pi-headAngle)*mechTrail)#
            turnAnglePnt2RearAxle = ((frontWheelRadius-turnAnglePntZ)**2+turnAnglePntX**2)**.5#
            
            
            ## front axle moves along the x-z plane with steer, sweeping through the torus
            ## steerAngle==pi at 180
            frontAxleZ = frontWheelRadius-((steerAngle/math.pi)*2)*(math.sin(math.pi/2-headAngle)*mechTrail)#
            frontAxleX = wheelBase-((steerAngle/math.pi)*2)*(math.sin(headAngle)*mechTrail)#
            ## angle between initial place, 0 steer, and current steer, frame not moved yet...
            ## relative to rear axle height; add any existing angle
            frontAxleAngle = math.atan((rearWheelRadius-frontAxleZ)/frontAxleX)+angleBtwnAxles#
            
            ## now, rotate the steer, calc angle due to steer without trail - holding frame fixed            
            ## max fWheel drop is prop. to steer Angle==90
            steerDropAngle = math.asin( (math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) /     \
                                       (wheelBase+math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) )#
            
            ## subtract the amount the wheel drops due to flop, from the amount 
            ## the frame rises due to the whole wheel sweeping about the center of the torus
            ## frontAxleAngle-steerDropAngle
            
            ## calc rearMassHeight
            ## so, hypotenous == rmass2Raxle
            newrmass2RaxleAngle = rmass2RaxleAngle+frontAxleAngle-steerDropAngle#
            newrearMassZ = rearWheelRadius+rmass2Raxle*math.sin(newrmass2RaxleAngle)#
                            
            #print '\t', round(newrearMassZ-rearMassZ, 2),
            #print '\t', round(rearPitchAngle/math.pi * 180.0, 2),
            #print ' ', round(frontAxle2Ground, 1),
            #print ' ', round(turnAnglePntZ, 1),
            print '\t', round(newrearMassZ-rearMassZ, 3),
            
if __name__ == '__main__':
    main()

JPEG image

#!/usr/bin/python
#Boa:PyApp:main

import math

modules ={}

def main():

    wheelBase = 72.5
    steerAngle = 20.  * math.pi / 180.0
    rearMassZ = 28.
    rearMassX = 67.
    frontWheelRadius = 33.
    rearWheelRadius = 33.


    ## calc angle from raxle to rmass, from horizontal
    rmass2RaxleAngle = math.asin((rearMassZ-rearWheelRadius)/rearMassX)
    ## calc dist from raxle to rmass
    rmass2Raxle = ((rearWheelRadius-rearMassZ)**2+rearMassX**2)**.5
    ## angle front to rear axle
    ## a wheelie is a pos+ angle
    angleBtwnAxles = math.atan((frontWheelRadius-rearWheelRadius)/wheelBase)
    
    print 'steerAngle', steerAngle/math.pi*180.0
    print 'pivot                           ground trail                      '
    print '     ', '\t10', '\t5', '\t0', '\t-5', '\t-10', '\t-15', '\t-20', 
'\t-25'
    for angle in range(40, 71, 1):        
        headAngle = angle  * math.pi / 180.0 # deg from horiz
        print '\n',int(round(headAngle*180/math.pi)), '   ',
        ## iterate over the ground trail
        for gtrail in range(10, -30, -5):
#            print "gtrail=", gtrail, "\n";
            mechTrail = frontWheelRadius*math.cos(headAngle) - 
gtrail*math.sin(headAngle) 
#            print "mechTrail=", mechTrail, "\n";
                
            ## rear wheel does not move or roll, frame rotates around the axle
            ## turnAnglePnt is imaginary point the axle rotates about with 
steer, 
            ## is fixed with respect to rear frame, and is the center of a torus
            ## swept by the wheel
            turnAnglePntZ = 
frontWheelRadius+(math.cos(math.pi-headAngle)*mechTrail)#
            turnAnglePntX = wheelBase-(math.sin(math.pi-headAngle)*mechTrail)#
            turnAnglePnt2RearAxle = 
((frontWheelRadius-turnAnglePntZ)**2+turnAnglePntX**2)**.5#
            
            
            ## front axle moves along the x-z plane with steer, sweeping 
through the torus
            ## steerAngle==pi at 180
            frontAxleZ = -math.cos(headAngle)*mechTrail*(1 - 
math.cos(steerAngle))
            frontAxleX = wheelBase - math.sin(headAngle)*mechTrail*(1 - 
math.cos(steerAngle))
#            print "frontAxleZ=", frontAxleZ, "\n"
#            print "frontAxleX=", frontAxleX, "\n"
            ## angle between initial place, 0 steer, and current steer, frame 
not moved yet...
            ## relative to rear axle height; add any existing angle

            ## now, rotate the steer, calc angle due to steer without trail - 
holding frame fixed            
            ## max fWheel drop is prop. to steer Angle==90
            SA=steerAngle
            HA=headAngle
            R=frontWheelRadius
            cosHA=math.cos(HA)
            sinHA=math.sin(HA)
            cosSA=math.cos(SA)
            sinSA=math.sin(HA)
            t1=(sinHA*sinHA + cosSA*(1-sinHA*sinHA))
            t2=(cosSA*cosHA*sinHA - cosHA*sinHA)
            if(t2 == 0):
                t=(-math.pi/2)
            else:
                t= math.atan(t1 / t2)
            Z=-cosHA*sinHA*R*math.cos(t) + cosSA*cosHA*sinHA*R*math.cos(t) + 
sinHA*sinHA*R*math.sin(t) + cosSA*(1-sinHA*sinHA)*R*math.sin(t)
#            print "Z=", Z, "\n";
            X=frontAxleX-cosHA*cosHA*R*math.cos(t) + 
cosSA*(1-cosHA*cosHA)*R*math.cos(t) -sinHA*cosHA*R*math.sin(t) + 
cosSA*sinHA*cosHA*R*math.sin(t)
#            print "X=", X, "\n";
            combinedDropAngle = math.atan((R+frontAxleZ+Z)/X)
            #steerDropAngle = math.asin( 
(math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) /   
(wheelBase+math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) )
            ## subtract the amount the wheel drops due to flop, from the amount 
            ## the frame rises due to the whole wheel sweeping about the center 
of the torus
            ## frontAxleAngle-steerDropAngle
            
            ## calc rearMassHeight
            ## so, hypotenous == rmass2Raxle
            newrmass2RaxleAngle = rmass2RaxleAngle-combinedDropAngle
            newrearMassZ = 
rearWheelRadius+rmass2Raxle*math.sin(newrmass2RaxleAngle)#
                            
            #print '\t', round(newrearMassZ-rearMassZ, 2),
            #print '\t', round(rearPitchAngle/math.pi * 180.0, 2),
            #print ' ', round(frontAxle2Ground, 1),
            #print ' ', round(turnAnglePntZ, 1),
            print '\t', round(newrearMassZ-rearMassZ, 3),
            
if __name__ == '__main__':
    main()

Other related posts: