eic-opticks/ana/cumtrapz.py

0001 #!/usr/bin/env python
0002 """
0003
0004 https://docs.sympy.org/latest/modules/functions/elementary.html#piecewise
0005
0006 Hmm seems that sympy doesnt like a mix of symbolic and floating point
0007 maybe better to use
0008
0009
0010
0011 """
0012 import logging
0013 log = logging.getLogger(__name__)
0014 import numpy as np
0015 import matplotlib.pyplot as plt
0016 import scipy.integrate as integrate
0017
0018 if __name__ == '__main__':
0019
0020     ri = np.array([
0021            [ 1.55 ,  1.478],
0022            [ 1.795,  1.48 ],
0023            [ 2.105,  1.484],
0024            [ 2.271,  1.486],
0025            [ 2.551,  1.492],
0026            [ 2.845,  1.496],
0027            [ 3.064,  1.499],
0028            [ 4.133,  1.526],
0029            [ 6.2  ,  1.619],
0030            [ 6.526,  1.618],
0031            [ 6.889,  1.527],
0032            [ 7.294,  1.554],
0033            [ 7.75 ,  1.793],
0034            [ 8.267,  1.783],
0035            [ 8.857,  1.664],
0036            [ 9.538,  1.554],
0037            [10.33 ,  1.454],
0038            [15.5  ,  1.454]
0039           ])
0040
0041
0042     en = ri[:,0]
0043     bi = 1.
0044     ct = bi/ri[:,1]
0045     s2 = (1.-ct)*(1.+ct)
0046
0047     s2i = integrate.cumtrapz( s2, en, initial=0. )  # initial keeps s2i shape same as s2 and en
0048     s2i_0 = integrate.cumtrapz( s2, en )
0049
0050     fig, ax = plt.subplots(figsize=[12.8, 7.2])
0051     ax.plot( en, s2i, label="s2i" )
0052     ax.scatter( en, s2i, label="s2i" )
0053
0054     ax.plot( en[:-1], s2i_0, label="s2i_0" )
0055
0056
0057     ax.legend()
0058     fig.show()
0059
0060
0061
0062