Let me take some time to try to explain what we did. I will have a look at
the Prop classes and try to see how things can be interfaced. 

I have another associate Patrick Hendriks who is by far the best of all of
us in OO and C++. He prepares a thesis on a general muon reconstruction
package. He has another year and a bit and now has to use it to study some
physics channel. For this package we wrote some classes to do tracking in
a magnetic field. To give you an idea, Patrick wrote more than 600
highly templated classes full of STL, the tracking in magnetic field are
only some 10 - 15 classes.

Some of our design ideas: a magnetic field is realy only needed by the
tracking and can therefor be hidden completely behind the track classes.
There is one method in the interface: to create it (only one, of course,
so it is a singleton). The real interface is in the TrackInMagneticField
class. That class can ask a service from the MagneticField class to return
a value for the field. You pass it a Point object and it returns <vector>
(Bx,By,Bz).

We realized the things you realy need in tracking is:
-1- extrapolate a track over a distance <step> in a magnetic field, or,
-2- extrapolate a track until it crosses a boundary (the next counter
plane, the next cylinder with silicon strips, the calorimeter, ...etc)
The second method obviously repeatedly uses the first one.
These are methods in the TrackInMagneticField class.

Realising that the thing you want to extrapolate to is polymorph (a plane,
a cylinder, a volume, ...etc) we wrote an abstract baseclass Interceptor
which implements the interface. The specific calculation to see if one has
crossed the boundary is done in the specialised classes VolumeInterceptor,
PlaneInterceptor, CylinderInterceptor. If you want yet something else you
have towrite your own YourInterceptor class, but the interface is the
same.

The MagneticField follows the same scheme: it's an abstract base class
with a unique public interface. There are several implementations like
ConstantMagneticField and FieldMap. If you want yet another field (like a
PolynomialMagneticField for example) you have to write it. This way the
interface is always the same but the calculation of the value of the field
to be returned could be completely different. This way it will be very
easy to switch from one sort of field to another without changing anything
else in the reconstruction code. This may be important for performance
reasons: looking up a value in a map takes much more time that calculating
from a polynomial.

The storage of the fieldmap is technically a bit complicated. It uses STL
<map>'s (you're now beginning to realise that I'm a great enthousiast of
STL, I supppose -:). In fact its a map_of_maps_of_maps: the first map
contains the x-values and each x-value points to another map which
contains all y-values etc. This makes the lookup of an (x,y,z) point very
fast. The last map contains z and (Bx,By,Bz). This bit is where the
reconstruction spends 70% of its time so we spend a lot of time to
optimise it for speed but making sure it cannot bumb. We needed maps to
allow for any grid in the map. The Saclay people provide denser grids
where the field changes steeply and a wide grid where the field is flat.
Such a non regular grid can not be accomodated with something simple as a
matrix. You need something fancier...and STL has it all!

There is of course also an interpolation machine. The field is rarely
asked for a value in a point which happens to be a grid point. For any
other point you have find the box of gridpoints surrounding it and do an
interpolation.

You imagine, the other classes which implement MagneticField
(ConstantMag.., PolynomialMagn.., DipoleMagn...etc) are much easier.

So what else is there to say? It's (mostly) written, there is some doc,
it's tested and it is used by Patrick now and Onne will probably use it
next. We're still trying to work out an algorith to calculate the stepsize
in the tracking. You want to make big steps where the field is smooth but
reduce the stepsize when you realise the field is steeper. We may deduce
it from the spacing of the grid. We haven't decide, we're still trying.

If you think, this could be usefull, I could present a bit more detail. I
was putting the code onto the web but couldn't finish it today. I'll
certainly put the doc out. Laurent Chevalier from Saclay participated in
the work and he's in the Saclay group which also joined D0 and took
responsibility for the magnetic field calculations for D0. So, there
should be no problem using their D0 fieldmap. Anyway, what we wrote is
very general and should fit any detector and any field and any geometry of
chambers....at least that's how we tried to write it -:)

I hope this helps,
cheers, 
kors

ps. 
I send a copy to John Womersley since this was discussed this morning
in the Algorithms meeting, but I couldn't see who was talking because the
camera was pointed at the screen and I didn't want to interupt.