Refactored obsolescent code

Fixes #19

Test case (same result from refactored version and the one in scipy):

 n    =           2
 mp   =           5
 x    =   1.2345600000000001     
 sj   =  0.76464742242968342       0.35211808607480077        9.1005002907604604E-002
 dj   = -0.35211808607480077       0.19421247220482474       0.13097451369045834     
 nm1  =           2
 sy   = -0.26725029005084994      -0.98112154280524222       -2.1168903174657765     
 dy   =  0.98112154280524222        1.3221767816269010        4.1629548507255132     
 nm2  =           2
 msta1=           6
 msta2=          16

This is the only one I did.

Here's another one (same result for both):

 n    =           4
 mp   =           6
 x    =   10.234560000000000     
 sj   =  -7.0754079663227554E-002   6.0471793911910328E-002   8.8479842347283441E-002  -1.7245782071985880E-002 -0.10027521747859354     
 dj   =  -6.0471793911910390E-002  -8.2571254785931469E-002   3.4536186807955657E-002   9.5220056708032075E-002   3.1742751621008022E-002
 nm1  =           4
 sy   =   6.7385044668486882E-002   7.7338148218054237E-002  -4.4715340777536736E-002  -9.9183416005597977E-002  -2.3121861020995477E-002
 dy   =  -7.7338148218054237E-002   5.2271908741186784E-002   9.0445308890580486E-002  -5.9512254640897572E-003  -8.7887443818715744E-002
 nm2  =           4
 msta1=          22
 msta2=          32

@jacobwilliams thanks for doing the tests. If you wouldn't mind doing a few more just to make sure, that would be great. If there is a bug introduced by this, we can always trace it to this PR and fix it.