shapes.py 124 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557155815591560156115621563156415651566156715681569157015711572157315741575157615771578157915801581158215831584158515861587158815891590159115921593159415951596159715981599160016011602160316041605160616071608160916101611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680168116821683168416851686168716881689169016911692169316941695169616971698169917001701170217031704170517061707170817091710171117121713171417151716171717181719172017211722172317241725172617271728172917301731173217331734173517361737173817391740174117421743174417451746174717481749175017511752175317541755175617571758175917601761176217631764176517661767176817691770177117721773177417751776177717781779178017811782178317841785178617871788178917901791179217931794179517961797179817991800180118021803180418051806180718081809181018111812181318141815181618171818181918201821182218231824182518261827182818291830183118321833183418351836183718381839184018411842184318441845184618471848184918501851185218531854185518561857185818591860186118621863186418651866186718681869187018711872187318741875187618771878187918801881188218831884188518861887188818891890189118921893189418951896189718981899190019011902190319041905190619071908190919101911191219131914191519161917191819191920192119221923192419251926192719281929193019311932193319341935193619371938193919401941194219431944194519461947194819491950195119521953195419551956195719581959196019611962196319641965196619671968196919701971197219731974197519761977197819791980198119821983198419851986198719881989199019911992199319941995199619971998199920002001200220032004200520062007200820092010201120122013201420152016201720182019202020212022202320242025202620272028202920302031203220332034203520362037203820392040204120422043204420452046204720482049205020512052205320542055205620572058205920602061206220632064206520662067206820692070207120722073207420752076207720782079208020812082208320842085208620872088208920902091209220932094209520962097209820992100210121022103210421052106210721082109211021112112211321142115211621172118211921202121212221232124212521262127212821292130213121322133213421352136213721382139214021412142214321442145214621472148214921502151215221532154215521562157215821592160216121622163216421652166216721682169217021712172217321742175217621772178217921802181218221832184218521862187218821892190219121922193219421952196219721982199220022012202220322042205220622072208220922102211221222132214221522162217221822192220222122222223222422252226222722282229223022312232223322342235223622372238223922402241224222432244224522462247224822492250225122522253225422552256225722582259226022612262226322642265226622672268226922702271227222732274227522762277227822792280228122822283228422852286228722882289229022912292229322942295229622972298229923002301230223032304230523062307230823092310231123122313231423152316231723182319232023212322232323242325232623272328232923302331233223332334233523362337233823392340234123422343234423452346234723482349235023512352235323542355235623572358235923602361236223632364236523662367236823692370237123722373237423752376237723782379238023812382238323842385238623872388238923902391239223932394239523962397239823992400240124022403240424052406240724082409241024112412241324142415241624172418241924202421242224232424242524262427242824292430243124322433243424352436243724382439244024412442244324442445244624472448244924502451245224532454245524562457245824592460246124622463246424652466246724682469247024712472247324742475247624772478247924802481248224832484248524862487248824892490249124922493249424952496249724982499250025012502250325042505250625072508250925102511251225132514251525162517251825192520252125222523252425252526252725282529253025312532253325342535253625372538253925402541254225432544254525462547254825492550255125522553255425552556255725582559256025612562256325642565256625672568256925702571257225732574257525762577257825792580258125822583258425852586258725882589259025912592259325942595259625972598259926002601260226032604260526062607260826092610261126122613261426152616261726182619262026212622262326242625262626272628262926302631263226332634263526362637263826392640
  1. from __future__ import division
  2. from __future__ import unicode_literals
  3. from __future__ import print_function
  4. from __future__ import absolute_import
  5. from future import standard_library
  6. standard_library.install_aliases()
  7. from builtins import input
  8. from builtins import zip
  9. from builtins import str
  10. from builtins import range
  11. from builtins import *
  12. from builtins import object
  13. from numpy import linspace, sin, cos, pi, array, asarray, ndarray, sqrt, abs
  14. import pprint, copy, glob, os
  15. from math import radians
  16. from .MatplotlibDraw import MatplotlibDraw
  17. drawing_tool = MatplotlibDraw()
  18. def point(x, y, check_inside=False):
  19. for obj, name in zip([x, y], ['x', 'y']):
  20. if isinstance(obj, (float,int)):
  21. pass
  22. elif isinstance(obj, ndarray):
  23. if obj.size == 1:
  24. pass
  25. else:
  26. raise TypeError('%s=%s of type %d has length=%d > 1' %
  27. (name, obj, type(obj), obj.size))
  28. else:
  29. raise TypeError('%s=%s is of wrong type %d' %
  30. (name, obj, type(obj)))
  31. if check_inside:
  32. ok, msg = drawing_tool.inside((x,y), exception=True)
  33. if not ok:
  34. print(msg)
  35. return array((x, y), dtype=float)
  36. def distance(p1, p2):
  37. p1 = arr2D(p1); p2 = arr2D(p2)
  38. d = p2 - p1
  39. return sqrt(d[0]**2 + d[1]**2)
  40. def unit_vec(x, y=None):
  41. """Return unit vector of the vector (x,y), or just x if x is a 2D point."""
  42. if isinstance(x, (float,int)) and isinstance(y, (float,int)):
  43. x = point(x, y)
  44. elif isinstance(x, (list,tuple,ndarray)) and y is None:
  45. return arr2D(x)/sqrt(x[0]**2 + x[1]**2)
  46. else:
  47. raise TypeError('x=%s is %s, must be float or ndarray 2D point' %
  48. (x, type(x)))
  49. def arr2D(x, check_inside=False):
  50. if isinstance(x, (tuple,list,ndarray)):
  51. if len(x) == 2:
  52. pass
  53. else:
  54. raise ValueError('x=%s has length %d, not 2' % (x, len(x)))
  55. else:
  56. raise TypeError('x=%s must be list/tuple/ndarray, not %s' %
  57. (x, type(x)))
  58. if check_inside:
  59. ok, msg = drawing_tool.inside(x, exception=True)
  60. if not ok:
  61. print(msg)
  62. return asarray(x, dtype=float)
  63. def _is_sequence(seq, length=None,
  64. can_be_None=False, error_message=True):
  65. if can_be_None:
  66. legal_types = (list,tuple,ndarray,None)
  67. else:
  68. legal_types = (list,tuple,ndarray)
  69. if isinstance(seq, legal_types):
  70. if length is not None:
  71. if length == len(seq):
  72. return True
  73. elif error_message:
  74. raise TypeError('%s is %s; must be %s of length %d' %
  75. (str(seq), type(seq),
  76. ', '.join([str(t) for t in legal_types]),
  77. len(seq)))
  78. else:
  79. return False
  80. else:
  81. return True
  82. elif error_message:
  83. raise TypeError('%s is %s, %s; must be %s' %
  84. (str(seq), seq.__class__.__name__, type(seq),
  85. ','.join([str(t)[5:-1] for t in legal_types])))
  86. else:
  87. return False
  88. def is_sequence(*sequences, **kwargs):
  89. length = kwargs.get('length', 2)
  90. can_be_None = kwargs.get('can_be_None', False)
  91. error_message = kwargs.get('error_message', True)
  92. check_inside = kwargs.get('check_inside', False)
  93. for x in sequences:
  94. _is_sequence(x, length=length, can_be_None=can_be_None,
  95. error_message=error_message)
  96. if check_inside:
  97. ok, msg = drawing_tool.inside(x, exception=True)
  98. if not ok:
  99. print(msg)
  100. def animate(fig, time_points, action, moviefiles=False,
  101. pause_per_frame=0.5, show_screen_graphics=True,
  102. title=None,
  103. **action_kwargs):
  104. if moviefiles:
  105. # Clean up old frame files
  106. framefilestem = 'tmp_frame_'
  107. framefiles = glob.glob('%s*.png' % framefilestem)
  108. for framefile in framefiles:
  109. os.remove(framefile)
  110. for n, t in enumerate(time_points):
  111. drawing_tool.erase()
  112. action(t, fig, **action_kwargs)
  113. #could demand returning fig, but in-place modifications
  114. #are done anyway
  115. #fig = action(t, fig)
  116. #if fig is None:
  117. # raise TypeError(
  118. # 'animate: action returns None, not fig\n'
  119. # '(a Shape object with the whole figure)')
  120. fig.draw()
  121. drawing_tool.display(title=title, show=show_screen_graphics)
  122. if moviefiles:
  123. drawing_tool.savefig('%s%04d.png' % (framefilestem, n))
  124. if moviefiles:
  125. return '%s%%04d.png' % framefilestem
  126. class Shape(object):
  127. """
  128. Superclass for drawing different geometric shapes.
  129. Subclasses define shapes, but drawing, rotation, translation,
  130. etc. are done in generic functions in this superclass.
  131. """
  132. def __init__(self):
  133. """
  134. Never to be called from subclasses.
  135. """
  136. raise NotImplementedError(
  137. 'class %s must implement __init__,\nwhich defines '
  138. 'self.shapes as a dict (or list) of Shape objects\n'
  139. 'Do not call Shape.__init__!' % \
  140. self.__class__.__name__)
  141. def set_name(self, name):
  142. self.name = name
  143. return self
  144. def get_name(self):
  145. return self.name if hasattr(self, 'name') else 'no_name'
  146. def __iter__(self):
  147. # We iterate over self.shapes many places, and will
  148. # get here if self.shapes is just a Shape object and
  149. # not the assumed dict/list.
  150. print('Warning: class %s does not define self.shapes\n'\
  151. 'as a dict of Shape objects')
  152. return [self] # Make the iteration work
  153. def copy(self):
  154. return copy.deepcopy(self)
  155. def __getitem__(self, name):
  156. """
  157. Allow indexing like::
  158. obj1['name1']['name2']
  159. all the way down to ``Curve`` or ``Point`` (``Text``)
  160. objects.
  161. """
  162. if hasattr(self, 'shapes'):
  163. if name in self.shapes:
  164. return self.shapes[name]
  165. else:
  166. for shape in self.shapes:
  167. if isinstance(self.shapes[shape], (Curve,Point)):
  168. # Indexing of Curve/Point/Text is not possible
  169. raise TypeError(
  170. 'Index "%s" (%s) is illegal' %
  171. (name, self.__class__.__name__))
  172. return self.shapes[shape][name]
  173. else:
  174. raise Exception('This is a bug in __getitem__')
  175. def __setitem__(self, name, value):
  176. """
  177. Allow assignment like::
  178. obj1['name1']['name2'] = value
  179. all the way down to ``Curve`` or ``Point`` (``Text``)
  180. objects.
  181. """
  182. if hasattr(self, 'shapes'):
  183. self.shapes[name] = value
  184. else:
  185. raise Exception('Cannot assign')
  186. def _for_all_shapes(self, func, *args, **kwargs):
  187. if not hasattr(self, 'shapes'):
  188. # When self.shapes is lacking, we either come to
  189. # a special implementation of func or we come here
  190. # because Shape.func is just inherited. This is
  191. # an error if the class is not Curve or Point
  192. if isinstance(self, (Curve, Point)):
  193. return # ok: no shapes and no func
  194. else:
  195. raise AttributeError('class %s has no shapes attribute!' %
  196. self.__class__.__name__)
  197. is_dict = True if isinstance(self.shapes, dict) else False
  198. for k, shape in enumerate(self.shapes):
  199. if is_dict:
  200. shape_name = shape
  201. shape = self.shapes[shape]
  202. else:
  203. shape_name = k # use index as name if list (not dict)
  204. if not isinstance(shape, Shape):
  205. if isinstance(shape, dict):
  206. raise TypeError(
  207. 'class %s has a self.shapes member "%s" that is just\n'
  208. 'a plain dictionary,\n%s\n'
  209. 'Did you mean to embed this dict in a Composition\n'
  210. 'object?' % (self.__class__.__name__, shape_name,
  211. str(shape)))
  212. elif isinstance(shape, (list,tuple)):
  213. raise TypeError(
  214. 'class %s has self.shapes member "%s" containing\n'
  215. 'a %s object %s,\n'
  216. 'Did you mean to embed this list in a Composition\n'
  217. 'object?' % (self.__class__.__name__, shape_name,
  218. type(shape), str(shape)))
  219. elif shape is None:
  220. raise TypeError(
  221. 'class %s has a self.shapes member "%s" that is None.\n'
  222. 'Some variable name is wrong, or some function\n'
  223. 'did not return the right object...' \
  224. % (self.__class__.__name__, shape_name))
  225. else:
  226. raise TypeError(
  227. 'class %s has a self.shapes member "%s" of %s which '
  228. 'is not a Shape object\n%s' %
  229. (self.__class__.__name__, shape_name, type(shape),
  230. pprint.pformat(self.shapes)))
  231. if isinstance(shape, Curve):
  232. shape.name = shape_name
  233. getattr(shape, func)(*args, **kwargs)
  234. def draw(self):
  235. self._for_all_shapes('draw')
  236. return self
  237. def draw_dimensions(self):
  238. if hasattr(self, 'dimensions'):
  239. for shape in self.dimensions:
  240. self.dimensions[shape].draw()
  241. return self
  242. else:
  243. #raise AttributeError('no self.dimensions dict for defining dimensions of class %s' % self.__classname__.__name__)
  244. return self
  245. def rotate(self, angle, center):
  246. is_sequence(center, length=2)
  247. self._for_all_shapes('rotate', angle, center)
  248. return self
  249. def translate(self, vec):
  250. is_sequence(vec, length=2)
  251. self._for_all_shapes('translate', vec)
  252. return self
  253. def scale(self, factor):
  254. self._for_all_shapes('scale', factor)
  255. return self
  256. def deform(self, displacement_function):
  257. self._for_all_shapes('deform', displacement_function)
  258. return self
  259. def minmax_coordinates(self, minmax=None):
  260. if minmax is None:
  261. minmax = {'xmin': 1E+20, 'xmax': -1E+20,
  262. 'ymin': 1E+20, 'ymax': -1E+20}
  263. self._for_all_shapes('minmax_coordinates', minmax)
  264. return minmax
  265. def recurse(self, name, indent=0):
  266. if not isinstance(self.shapes, dict):
  267. raise TypeError('recurse works only with dict self.shape, not %s' %
  268. type(self.shapes))
  269. space = ' '*indent
  270. print(space, '%s: %s.shapes has entries' % \
  271. (self.__class__.__name__, name), \
  272. str(list(self.shapes.keys()))[1:-1])
  273. for shape in self.shapes:
  274. print(space, end=' ')
  275. print('call %s.shapes["%s"].recurse("%s", %d)' % \
  276. (name, shape, shape, indent+2))
  277. self.shapes[shape].recurse(shape, indent+2)
  278. def graphviz_dot(self, name, classname=True):
  279. if not isinstance(self.shapes, dict):
  280. raise TypeError('recurse works only with dict self.shape, not %s' %
  281. type(self.shapes))
  282. dotfile = name + '.dot'
  283. pngfile = name + '.png'
  284. if classname:
  285. name = r"%s:\n%s" % (self.__class__.__name__, name)
  286. couplings = self._object_couplings(name, classname=classname)
  287. # Insert counter for similar names
  288. from collections import defaultdict
  289. count = defaultdict(lambda: 0)
  290. couplings2 = []
  291. for i in range(len(couplings)):
  292. parent, child = couplings[i]
  293. count[child] += 1
  294. parent += ' (%d)' % count[parent]
  295. child += ' (%d)' % count[child]
  296. couplings2.append((parent, child))
  297. print('graphviz', couplings, count)
  298. # Remove counter for names there are only one of
  299. for i in range(len(couplings)):
  300. parent2, child2 = couplings2[i]
  301. parent, child = couplings[i]
  302. if count[parent] > 1:
  303. parent = parent2
  304. if count[child] > 1:
  305. child = child2
  306. couplings[i] = (parent, child)
  307. print(couplings)
  308. f = open(dotfile, 'w')
  309. f.write('digraph G {\n')
  310. for parent, child in couplings:
  311. f.write('"%s" -> "%s";\n' % (parent, child))
  312. f.write('}\n')
  313. f.close()
  314. print('Run dot -Tpng -o %s %s' % (pngfile, dotfile))
  315. def _object_couplings(self, parent, couplings=[], classname=True):
  316. """Find all couplings of parent and child objects in a figure."""
  317. for shape in self.shapes:
  318. if classname:
  319. childname = r"%s:\n%s" % \
  320. (self.shapes[shape].__class__.__name__, shape)
  321. else:
  322. childname = shape
  323. couplings.append((parent, childname))
  324. self.shapes[shape]._object_couplings(childname, couplings,
  325. classname)
  326. return couplings
  327. def set_linestyle(self, style):
  328. styles = ('solid', 'dashed', 'dashdot', 'dotted')
  329. if style not in styles:
  330. raise ValueError('%s: style=%s must be in %s' %
  331. (self.__class__.__name__ + '.set_linestyle:',
  332. style, str(styles)))
  333. self._for_all_shapes('set_linestyle', style)
  334. return self
  335. def set_linewidth(self, width):
  336. if not isinstance(width, int) and width >= 0:
  337. raise ValueError('%s: width=%s must be positive integer' %
  338. (self.__class__.__name__ + '.set_linewidth:',
  339. width))
  340. self._for_all_shapes('set_linewidth', width)
  341. return self
  342. def set_linecolor(self, color):
  343. if color in drawing_tool.line_colors:
  344. color = drawing_tool.line_colors[color]
  345. elif color in list(drawing_tool.line_colors.values()):
  346. pass # color is ok
  347. else:
  348. raise ValueError('%s: invalid color "%s", must be in %s' %
  349. (self.__class__.__name__ + '.set_linecolor:',
  350. color, list(drawing_tool.line_colors.keys())))
  351. self._for_all_shapes('set_linecolor', color)
  352. return self
  353. def set_arrow(self, style):
  354. styles = ('->', '<-', '<->')
  355. if not style in styles:
  356. raise ValueError('%s: style=%s must be in %s' %
  357. (self.__class__.__name__ + '.set_arrow:',
  358. style, styles))
  359. self._for_all_shapes('set_arrow', style)
  360. return self
  361. def set_filled_curves(self, color='', pattern=''):
  362. if color in drawing_tool.line_colors:
  363. color = drawing_tool.line_colors[color]
  364. elif color in list(drawing_tool.line_colors.values()):
  365. pass # color is ok
  366. else:
  367. raise ValueError('%s: invalid color "%s", must be in %s' %
  368. (self.__class__.__name__ + '.set_filled_curves:',
  369. color, list(drawing_tool.line_colors.keys())))
  370. self._for_all_shapes('set_filled_curves', color, pattern)
  371. return self
  372. def set_shadow(self, pixel_displacement=3):
  373. self._for_all_shapes('set_shadow', pixel_displacement)
  374. return self
  375. def show_hierarchy(self, indent=0, format='std'):
  376. """Recursive pretty print of hierarchy of objects."""
  377. if not isinstance(self.shapes, dict):
  378. print('cannot print hierarchy when %s.shapes is not a dict' % \
  379. self.__class__.__name__)
  380. s = ''
  381. if format == 'dict':
  382. s += '{'
  383. for shape in self.shapes:
  384. if format == 'dict':
  385. shape_str = repr(shape) + ':'
  386. elif format == 'plain':
  387. shape_str = shape
  388. else:
  389. shape_str = shape + ':'
  390. if format == 'dict' or format == 'plain':
  391. class_str = ''
  392. else:
  393. class_str = ' (%s)' % \
  394. self.shapes[shape].__class__.__name__
  395. s += '\n%s%s%s %s,' % (
  396. ' '*indent,
  397. shape_str,
  398. class_str,
  399. self.shapes[shape].show_hierarchy(indent+4, format))
  400. if format == 'dict':
  401. s += '}'
  402. return s
  403. def __str__(self):
  404. """Display hierarchy with minimum information (just object names)."""
  405. return self.show_hierarchy(format='plain')
  406. def __repr__(self):
  407. """Display hierarchy as a dictionary."""
  408. return self.show_hierarchy(format='dict')
  409. #return pprint.pformat(self.shapes)
  410. class Curve(Shape):
  411. """General curve as a sequence of (x,y) coordintes."""
  412. def __init__(self, x, y):
  413. """
  414. `x`, `y`: arrays holding the coordinates of the curve.
  415. """
  416. self.x = asarray(x, dtype=float)
  417. self.y = asarray(y, dtype=float)
  418. #self.shapes must not be defined in this class
  419. #as self.shapes holds children objects:
  420. #Curve has no children (end leaf of self.shapes tree)
  421. self.linestyle = None
  422. self.linewidth = None
  423. self.linecolor = None
  424. self.fillcolor = None
  425. self.fillpattern = None
  426. self.arrow = None
  427. self.shadow = False
  428. self.name = None # name of object that this Curve represents
  429. def inside_plot_area(self, verbose=True):
  430. """Check that all coordinates are within drawing_tool's area."""
  431. xmin, xmax = self.x.min(), self.x.max()
  432. ymin, ymax = self.y.min(), self.y.max()
  433. t = drawing_tool
  434. inside = True
  435. if not hasattr(t, 'xmin'):
  436. return None # drawing area is not defined
  437. if xmin < t.xmin:
  438. inside = False
  439. if verbose:
  440. print('x_min=%g < plot area x_min=%g' % (xmin, t.xmin))
  441. if xmax > t.xmax:
  442. inside = False
  443. if verbose:
  444. print('x_max=%g > plot area x_max=%g' % (xmax, t.xmax))
  445. if ymin < t.ymin:
  446. inside = False
  447. if verbose:
  448. print('y_min=%g < plot area y_min=%g' % (ymin, t.ymin))
  449. if ymax > t.ymax:
  450. inside = False
  451. if verbose:
  452. print('y_max=%g > plot area y_max=%g' % (ymax, t.ymax))
  453. return inside
  454. def draw(self):
  455. """
  456. Send the curve to the plotting engine. That is, convert
  457. coordinate information in self.x and self.y, together
  458. with optional settings of linestyles, etc., to
  459. plotting commands for the chosen engine.
  460. """
  461. self.inside_plot_area()
  462. drawing_tool.plot_curve(
  463. self.x, self.y,
  464. self.linestyle, self.linewidth, self.linecolor,
  465. self.arrow, self.fillcolor, self.fillpattern,
  466. self.shadow, self.name)
  467. def rotate(self, angle, center):
  468. """
  469. Rotate all coordinates: `angle` is measured in degrees and
  470. (`x`,`y`) is the "origin" of the rotation.
  471. """
  472. angle = radians(angle)
  473. x, y = center
  474. c = cos(angle); s = sin(angle)
  475. xnew = x + (self.x - x)*c - (self.y - y)*s
  476. ynew = y + (self.x - x)*s + (self.y - y)*c
  477. self.x = xnew
  478. self.y = ynew
  479. return self
  480. def scale(self, factor):
  481. """Scale all coordinates by `factor`: ``x = factor*x``, etc."""
  482. self.x = factor*self.x
  483. self.y = factor*self.y
  484. return self
  485. def translate(self, vec):
  486. """Translate all coordinates by a vector `vec`."""
  487. self.x += vec[0]
  488. self.y += vec[1]
  489. return self
  490. def deform(self, displacement_function):
  491. """Displace all coordinates according to displacement_function(x,y)."""
  492. for i in range(len(self.x)):
  493. self.x[i], self.y[i] = displacement_function(self.x[i], self.y[i])
  494. return self
  495. def minmax_coordinates(self, minmax=None):
  496. if minmax is None:
  497. minmax = {'xmin': [], 'xmax': [], 'ymin': [], 'ymax': []}
  498. minmax['xmin'] = min(self.x.min(), minmax['xmin'])
  499. minmax['xmax'] = max(self.x.max(), minmax['xmax'])
  500. minmax['ymin'] = min(self.y.min(), minmax['ymin'])
  501. minmax['ymax'] = max(self.y.max(), minmax['ymax'])
  502. return minmax
  503. def recurse(self, name, indent=0):
  504. space = ' '*indent
  505. print(space, 'reached "bottom" object %s' % \
  506. self.__class__.__name__)
  507. def _object_couplings(self, parent, couplings=[], classname=True):
  508. return
  509. def set_linecolor(self, color):
  510. self.linecolor = color
  511. return self
  512. def set_linewidth(self, width):
  513. self.linewidth = width
  514. return self
  515. def set_linestyle(self, style):
  516. self.linestyle = style
  517. return self
  518. def set_arrow(self, style=None):
  519. self.arrow = style
  520. return self
  521. def set_filled_curves(self, color='', pattern=''):
  522. self.fillcolor = color
  523. self.fillpattern = pattern
  524. return self
  525. def set_shadow(self, pixel_displacement=3):
  526. self.shadow = pixel_displacement
  527. return self
  528. def show_hierarchy(self, indent=0, format='std'):
  529. if format == 'dict':
  530. return '"%s"' % str(self)
  531. elif format == 'plain':
  532. return ''
  533. else:
  534. return str(self)
  535. def __str__(self):
  536. """Compact pretty print of a Curve object."""
  537. s = '%d (x,y) coords' % self.x.size
  538. inside = self.inside_plot_area(verbose=False)
  539. if inside is None:
  540. pass # no info about the plotting area
  541. elif not inside:
  542. s += ', some coordinates are outside plotting area!\n'
  543. props = ('linecolor', 'linewidth', 'linestyle', 'arrow',
  544. 'fillcolor', 'fillpattern')
  545. for prop in props:
  546. value = getattr(self, prop)
  547. if value is not None:
  548. s += ' %s=%s' % (prop, repr(value))
  549. return s
  550. def __repr__(self):
  551. return str(self)
  552. class Spline(Shape):
  553. # Note: UnivariateSpline interpolation may not work if
  554. # the x[i] points are far from uniformly spaced
  555. def __init__(self, x, y, degree=3, resolution=501):
  556. from scipy.interpolate import UnivariateSpline
  557. self.smooth = UnivariateSpline(x, y, s=0, k=degree)
  558. self.xcoor = linspace(x[0], x[-1], resolution)
  559. ycoor = self.smooth(self.xcoor)
  560. self.shapes = {'smooth': Curve(self.xcoor, ycoor)}
  561. def geometric_features(self):
  562. s = self.shapes['smooth']
  563. return {'start': point(s.x[0], s.y[0]),
  564. 'end': point(s.x[-1], s.y[-1]),
  565. 'interval': [s.x[0], s.x[-1]]}
  566. def __call__(self, x):
  567. return self.smooth(x)
  568. # Can easily find the derivative and the integral as
  569. # self.smooth.derivative(n=1) and self.smooth.antiderivative()
  570. class SketchyFunc1(Spline):
  571. """
  572. A typical function curve used to illustrate an "arbitrary" function.
  573. """
  574. domain = [1, 6]
  575. def __init__(self, name=None, name_pos='start',
  576. xmin=1, xmax=6, ymin=2.4, ymax=5):
  577. x = array([1, 2, 3, 4, 5, 6])
  578. y = array([5, 3.5, 3.8, 3, 2.5, 2.4])
  579. # Scale x and y
  580. x = xmin - x.min() + x*(xmax - xmin)/(x.max()-x.min())
  581. y = ymin - y.min() + y*(ymax - ymin)/(y.max()-y.min())
  582. Spline.__init__(self, x, y)
  583. self.shapes['smooth'].set_linecolor('black')
  584. if name is not None:
  585. self.shapes['name'] = Text(name, self.geometric_features()[name_pos] + point(0,0.1))
  586. class SketchyFunc3(Spline):
  587. """
  588. A typical function curve used to illustrate an "arbitrary" function.
  589. """
  590. domain = [0, 6]
  591. def __init__(self, name=None, name_pos='start',
  592. xmin=0, xmax=6, ymin=0.5, ymax=3.8):
  593. x = array([0, 2, 3, 4, 5, 6])
  594. #y = array([2, 3.5, 3.8, 2, 2.5, 2.6])
  595. y = array([0.5, 3.5, 3.8, 2, 2.5, 3.5])
  596. # Scale x and y
  597. x = xmin - x.min() + x*(xmax - xmin)/(x.max()-x.min())
  598. y = ymin - y.min() + y*(ymax - ymin)/(y.max()-y.min())
  599. Spline.__init__(self, x, y)
  600. self.shapes['smooth'].set_linecolor('black')
  601. if name is not None:
  602. self.shapes['name'] = Text(name, self.geometric_features()[name_pos] + point(0,0.1))
  603. class SketchyFunc4(Spline):
  604. """
  605. A typical function curve used to illustrate an "arbitrary" function.
  606. Can be a companion function to SketchyFunc3.
  607. """
  608. domain = [1, 6]
  609. def __init__(self, name=None, name_pos='start',
  610. xmin=0, xmax=6, ymin=0.5, ymax=1.8):
  611. x = array([0, 2, 3, 4, 5, 6])
  612. y = array([1.5, 1.3, 0.7, 0.5, 0.6, 0.8])
  613. # Scale x and y
  614. x = xmin - x.min() + x*(xmax - xmin)/(x.max()-x.min())
  615. y = ymin - y.min() + y*(ymax - ymin)/(y.max()-y.min())
  616. Spline.__init__(self, x, y)
  617. self.shapes['smooth'].set_linecolor('black')
  618. if name is not None:
  619. self.shapes['name'] = Text(name, self.geometric_features()[name_pos] + point(0,0.1))
  620. class SketchyFunc2(Shape):
  621. """
  622. A typical function curve used to illustrate an "arbitrary" function.
  623. """
  624. domain = [0, 2.25]
  625. def __init__(self, name=None, name_pos='end',
  626. xmin=0, xmax=2.25, ymin=0.046679703125, ymax=1.259375):
  627. a = 0; b = 2.25
  628. resolution = 100
  629. x = linspace(a, b, resolution+1)
  630. f = self # for calling __call__
  631. y = f(x)
  632. # Scale x and y
  633. x = xmin - x.min() + x*(xmax - xmin)/(x.max()-x.min())
  634. y = ymin - y.min() + y*(ymax - ymin)/(y.max()-y.min())
  635. self.shapes = {'smooth': Curve(x, y)}
  636. self.shapes['smooth'].set_linecolor('black')
  637. pos = point(a, f(a)) if name_pos == 'start' else point(b, f(b))
  638. if name is not None:
  639. self.shapes['name'] = Text(name, pos + point(0,0.1))
  640. def __call__(self, x):
  641. return 0.5+x*(2-x)*(0.9-x) # on [0, 2.25]
  642. class Point(Shape):
  643. """A point (x,y) which can be rotated, translated, and scaled."""
  644. def __init__(self, x, y):
  645. self.x, self.y = x, y
  646. #self.shapes is not needed in this class
  647. def __add__(self, other):
  648. if isinstance(other, (list,tuple)):
  649. other = Point(other)
  650. return Point(self.x+other.x, self.y+other.y)
  651. # class Point is an abstract class - only subclasses are useful
  652. # and must implement draw
  653. def draw(self):
  654. raise NotImplementedError(
  655. 'class %s must implement the draw method' %
  656. self.__class__.__name__)
  657. def rotate(self, angle, center):
  658. """Rotate point an `angle` (in degrees) around (`x`,`y`)."""
  659. angle = angle*pi/180
  660. x, y = center
  661. c = cos(angle); s = sin(angle)
  662. xnew = x + (self.x - x)*c - (self.y - y)*s
  663. ynew = y + (self.x - x)*s + (self.y - y)*c
  664. self.x = xnew
  665. self.y = ynew
  666. return self
  667. def scale(self, factor):
  668. """Scale point coordinates by `factor`: ``x = factor*x``, etc."""
  669. self.x = factor*self.x
  670. self.y = factor*self.y
  671. return self
  672. def translate(self, vec):
  673. """Translate point by a vector `vec`."""
  674. self.x += vec[0]
  675. self.y += vec[1]
  676. return self
  677. def deform(self, displacement_function):
  678. """Displace coordinates according to displacement_function(x,y)."""
  679. for i in range(len(self.x)):
  680. self.x, self.y = displacement_function(self.x, self.y)
  681. return self
  682. def minmax_coordinates(self, minmax=None):
  683. if minmax is None:
  684. minmax = {'xmin': [], 'xmax': [], 'ymin': [], 'ymax': []}
  685. minmax['xmin'] = min(self.x, minmax['xmin'])
  686. minmax['xmax'] = max(self.x, minmax['xmax'])
  687. minmax['ymin'] = min(self.y, minmax['ymin'])
  688. minmax['ymax'] = max(self.y, minmax['ymax'])
  689. return minmax
  690. def recurse(self, name, indent=0):
  691. space = ' '*indent
  692. print(space, 'reached "bottom" object %s' % \
  693. self.__class__.__name__)
  694. def _object_couplings(self, parent, couplings=[], classname=True):
  695. return
  696. # No need for set_linecolor etc since self._for_all_shapes, which
  697. # is always called for these functions, makes a test and stops
  698. # calls if self.shapes is missing and the object is Point or Curve
  699. def show_hierarchy(self, indent=0, format='std'):
  700. s = '%s at (%g,%g)' % (self.__class__.__name__, self.x, self.y)
  701. if format == 'dict':
  702. return '"%s"' % s
  703. elif format == 'plain':
  704. return ''
  705. else:
  706. return s
  707. # no need to store input data as they are invalid after rotations etc.
  708. class Rectangle(Shape):
  709. """
  710. Rectangle specified by the point `lower_left_corner`, `width`,
  711. and `height`.
  712. """
  713. def __init__(self, lower_left_corner, width, height):
  714. is_sequence(lower_left_corner)
  715. p = arr2D(lower_left_corner) # short form
  716. x = [p[0], p[0] + width,
  717. p[0] + width, p[0], p[0]]
  718. y = [p[1], p[1], p[1] + height,
  719. p[1] + height, p[1]]
  720. self.shapes = {'rectangle': Curve(x,y)}
  721. # Dimensions
  722. dims = {
  723. 'width': Distance_wText(p + point(0, -height/5.),
  724. p + point(width, -height/5.),
  725. 'width'),
  726. 'height': Distance_wText(p + point(width + width/5., 0),
  727. p + point(width + width/5., height),
  728. 'height'),
  729. 'lower_left_corner': Text_wArrow('lower_left_corner',
  730. p - point(width/5., height/5.), p)
  731. }
  732. self.dimensions = dims
  733. def geometric_features(self):
  734. """
  735. Return dictionary with
  736. ==================== =============================================
  737. Attribute Description
  738. ==================== =============================================
  739. lower_left Lower left corner point.
  740. upper_left Upper left corner point.
  741. lower_right Lower right corner point.
  742. upper_right Upper right corner point.
  743. lower_mid Middle point on lower side.
  744. upper_mid Middle point on upper side.
  745. center Center point
  746. ==================== =============================================
  747. """
  748. r = self.shapes['rectangle']
  749. d = {'lower_left': point(r.x[0], r.y[0]),
  750. 'lower_right': point(r.x[1], r.y[1]),
  751. 'upper_right': point(r.x[2], r.y[2]),
  752. 'upper_left': point(r.x[3], r.y[3])}
  753. d['lower_mid'] = 0.5*(d['lower_left'] + d['lower_right'])
  754. d['upper_mid'] = 0.5*(d['upper_left'] + d['upper_right'])
  755. d['left_mid'] = 0.5*(d['lower_left'] + d['upper_left'])
  756. d['right_mid'] = 0.5*(d['lower_right'] + d['upper_right'])
  757. d['center'] = point(d['lower_mid'][0], d['left_mid'][1])
  758. return d
  759. class Triangle(Shape):
  760. """
  761. Triangle defined by its three vertices p1, p2, and p3.
  762. Recorded geometric features:
  763. ==================== =============================================
  764. Attribute Description
  765. ==================== =============================================
  766. p1, p2, p3 Corners as given to the constructor.
  767. ==================== =============================================
  768. """
  769. def __init__(self, p1, p2, p3):
  770. is_sequence(p1, p2, p3)
  771. x = [p1[0], p2[0], p3[0], p1[0]]
  772. y = [p1[1], p2[1], p3[1], p1[1]]
  773. self.shapes = {'triangle': Curve(x,y)}
  774. # Dimensions
  775. self.dimensions = {'p1': Text('p1', p1),
  776. 'p2': Text('p2', p2),
  777. 'p3': Text('p3', p3)}
  778. def geometric_features(self):
  779. t = self.shapes['triangle']
  780. return {'p1': point(t.x[0], t.y[0]),
  781. 'p2': point(t.x[1], t.y[1]),
  782. 'p3': point(t.x[2], t.y[2])}
  783. class Line(Shape):
  784. def __init__(self, start, end):
  785. is_sequence(start, end, length=2)
  786. x = [start[0], end[0]]
  787. y = [start[1], end[1]]
  788. self.shapes = {'line': Curve(x, y)}
  789. def geometric_features(self):
  790. line = self.shapes['line']
  791. return {'start': point(line.x[0], line.y[0]),
  792. 'end': point(line.x[1], line.y[1]),}
  793. def compute_formulas(self):
  794. x, y = self.shapes['line'].x, self.shapes['line'].y
  795. # Define equations for line:
  796. # y = a*x + b, x = c*y + d
  797. try:
  798. self.a = (y[1] - y[0])/(x[1] - x[0])
  799. self.b = y[0] - self.a*x[0]
  800. except ZeroDivisionError:
  801. # Vertical line, y is not a function of x
  802. self.a = None
  803. self.b = None
  804. try:
  805. if self.a is None:
  806. self.c = 0
  807. else:
  808. self.c = 1/float(self.a)
  809. if self.b is None:
  810. self.d = x[1]
  811. except ZeroDivisionError:
  812. # Horizontal line, x is not a function of y
  813. self.c = None
  814. self.d = None
  815. def compute_formulas(self):
  816. x, y = self.shapes['line'].x, self.shapes['line'].y
  817. tol = 1E-14
  818. # Define equations for line:
  819. # y = a*x + b, x = c*y + d
  820. if abs(x[1] - x[0]) > tol:
  821. self.a = (y[1] - y[0])/(x[1] - x[0])
  822. self.b = y[0] - self.a*x[0]
  823. else:
  824. # Vertical line, y is not a function of x
  825. self.a = None
  826. self.b = None
  827. if self.a is None:
  828. self.c = 0
  829. elif abs(self.a) > tol:
  830. self.c = 1/float(self.a)
  831. self.d = x[1]
  832. else: # self.a is 0
  833. # Horizontal line, x is not a function of y
  834. self.c = None
  835. self.d = None
  836. def __call__(self, x=None, y=None):
  837. """Given x, return y on the line, or given y, return x."""
  838. self.compute_formulas()
  839. if x is not None and self.a is not None:
  840. return self.a*x + self.b
  841. elif y is not None and self.c is not None:
  842. return self.c*y + self.d
  843. else:
  844. raise ValueError(
  845. 'Line.__call__(x=%s, y=%s) not meaningful' % \
  846. (x, y))
  847. def new_interval(self, x=None, y=None):
  848. """Redefine current Line to cover interval in x or y."""
  849. if x is not None:
  850. is_sequence(x, length=2)
  851. xL, xR = x
  852. new_line = Line((xL, self(x=xL)), (xR, self(x=xR)))
  853. elif y is not None:
  854. is_sequence(y, length=2)
  855. yL, yR = y
  856. new_line = Line((xL, self(y=xL)), (xR, self(y=xR)))
  857. self.shapes['line'] = new_line['line']
  858. return self
  859. # First implementation of class Circle
  860. class Circle(Shape):
  861. def __init__(self, center, radius, resolution=180):
  862. self.center, self.radius = center, radius
  863. self.resolution = resolution
  864. t = linspace(0, 2*pi, resolution+1)
  865. x0 = center[0]; y0 = center[1]
  866. R = radius
  867. x = x0 + R*cos(t)
  868. y = y0 + R*sin(t)
  869. self.shapes = {'circle': Curve(x, y)}
  870. def __call__(self, theta):
  871. """
  872. Return (x, y) point corresponding to angle theta.
  873. Not valid after a translation, rotation, or scaling.
  874. """
  875. return self.center[0] + self.radius*cos(theta), \
  876. self.center[1] + self.radius*sin(theta)
  877. class Arc(Shape):
  878. def __init__(self, center, radius,
  879. start_angle, arc_angle,
  880. resolution=180):
  881. is_sequence(center)
  882. # Must record some parameters for __call__
  883. self.center = arr2D(center)
  884. self.radius = radius
  885. self.start_angle = radians(start_angle)
  886. self.arc_angle = radians(arc_angle)
  887. t = linspace(self.start_angle,
  888. self.start_angle + self.arc_angle,
  889. resolution+1)
  890. x0 = center[0]; y0 = center[1]
  891. R = radius
  892. x = x0 + R*cos(t)
  893. y = y0 + R*sin(t)
  894. self.shapes = {'arc': Curve(x, y)}
  895. # Cannot set dimensions (Arc_wText recurses into this
  896. # constructor forever). Set in test_Arc instead.
  897. def geometric_features(self):
  898. a = self.shapes['arc']
  899. m = len(a.x)//2 # mid point in array
  900. d = {'start': point(a.x[0], a.y[0]),
  901. 'end': point(a.x[-1], a.y[-1]),
  902. 'mid': point(a.x[m], a.y[m])}
  903. return d
  904. def __call__(self, theta):
  905. """
  906. Return (x,y) point at start_angle + theta.
  907. Not valid after translation, rotation, or scaling.
  908. """
  909. theta = radians(theta)
  910. t = self.start_angle + theta
  911. x0 = self.center[0]
  912. y0 = self.center[1]
  913. R = self.radius
  914. x = x0 + R*cos(t)
  915. y = y0 + R*sin(t)
  916. return (x, y)
  917. # Alternative for small arcs: Parabola
  918. class Parabola(Shape):
  919. def __init__(self, start, mid, stop, resolution=21):
  920. self.p1, self.p2, self.p3 = start, mid, stop
  921. # y as function of x? (no point on line x=const?)
  922. tol = 1E-14
  923. if abs(self.p1[0] - self.p2[0]) > 1E-14 and \
  924. abs(self.p2[0] - self.p3[0]) > 1E-14 and \
  925. abs(self.p3[0] - self.p1[0]) > 1E-14:
  926. self.y_of_x = True
  927. else:
  928. self.y_of_x = False
  929. # x as function of y? (no point on line y=const?)
  930. tol = 1E-14
  931. if abs(self.p1[1] - self.p2[1]) > 1E-14 and \
  932. abs(self.p2[1] - self.p3[1]) > 1E-14 and \
  933. abs(self.p3[1] - self.p1[1]) > 1E-14:
  934. self.x_of_y = True
  935. else:
  936. self.x_of_y = False
  937. if self.y_of_x:
  938. x = linspace(start[0], end[0], resolution)
  939. y = self(x=x)
  940. elif self.x_of_y:
  941. y = linspace(start[1], end[1], resolution)
  942. x = self(y=y)
  943. else:
  944. raise ValueError(
  945. 'Parabola: two or more points lie on x=const '
  946. 'or y=const - not allowed')
  947. self.shapes = {'parabola': Curve(x, y)}
  948. def __call__(self, x=None, y=None):
  949. if x is not None and self.y_of_x:
  950. return self._L2x(self.p1, self.p2)*self.p3[1] + \
  951. self._L2x(self.p2, self.p3)*self.p1[1] + \
  952. self._L2x(self.p3, self.p1)*self.p2[1]
  953. elif y is not None and self.x_of_y:
  954. return self._L2y(self.p1, self.p2)*self.p3[0] + \
  955. self._L2y(self.p2, self.p3)*self.p1[0] + \
  956. self._L2y(self.p3, self.p1)*self.p2[0]
  957. else:
  958. raise ValueError(
  959. 'Parabola.__call__(x=%s, y=%s) not meaningful' % \
  960. (x, y))
  961. def _L2x(self, x, pi, pj, pk):
  962. return (x - pi[0])*(x - pj[0])/((pk[0] - pi[0])*(pk[0] - pj[0]))
  963. def _L2y(self, y, pi, pj, pk):
  964. return (y - pi[1])*(y - pj[1])/((pk[1] - pi[1])*(pk[1] - pj[1]))
  965. class Circle(Arc):
  966. def __init__(self, center, radius, resolution=180):
  967. Arc.__init__(self, center, radius, 0, 360, resolution)
  968. class Wall(Shape):
  969. def __init__(self, x, y, thickness, pattern='/', transparent=False):
  970. is_sequence(x, y, length=len(x))
  971. if isinstance(x[0], (tuple,list,ndarray)):
  972. # x is list of curves
  973. x1 = concatenate(x)
  974. else:
  975. x1 = asarray(x, float)
  976. if isinstance(y[0], (tuple,list,ndarray)):
  977. # x is list of curves
  978. y1 = concatenate(y)
  979. else:
  980. y1 = asarray(y, float)
  981. self.x1 = x1; self.y1 = y1
  982. # Displaced curve (according to thickness)
  983. x2 = x1
  984. y2 = y1 + thickness
  985. # Combine x1,y1 with x2,y2 reversed
  986. from numpy import concatenate
  987. x = concatenate((x1, x2[-1::-1]))
  988. y = concatenate((y1, y2[-1::-1]))
  989. wall = Curve(x, y)
  990. wall.set_filled_curves(color='white', pattern=pattern)
  991. x = [x1[-1]] + x2[-1::-1].tolist() + [x1[0]]
  992. y = [y1[-1]] + y2[-1::-1].tolist() + [y1[0]]
  993. self.shapes = {'wall': wall}
  994. from collections import OrderedDict
  995. self.shapes = OrderedDict()
  996. self.shapes['wall'] = wall
  997. if transparent:
  998. white_eraser = Curve(x, y)
  999. white_eraser.set_linecolor('white')
  1000. self.shapes['eraser'] = white_eraser
  1001. def geometric_features(self):
  1002. d = {'start': point(self.x1[0], self.y1[0]),
  1003. 'end': point(self.x1[-1], self.y1[-1])}
  1004. return d
  1005. class Wall2(Shape):
  1006. def __init__(self, x, y, thickness, pattern='/'):
  1007. is_sequence(x, y, length=len(x))
  1008. if isinstance(x[0], (tuple,list,ndarray)):
  1009. # x is list of curves
  1010. x1 = concatenate(x)
  1011. else:
  1012. x1 = asarray(x, float)
  1013. if isinstance(y[0], (tuple,list,ndarray)):
  1014. # x is list of curves
  1015. y1 = concatenate(y)
  1016. else:
  1017. y1 = asarray(y, float)
  1018. self.x1 = x1; self.y1 = y1
  1019. # Displaced curve (according to thickness)
  1020. x2 = x1.copy()
  1021. y2 = y1.copy()
  1022. def displace(idx, idx_m, idx_p):
  1023. # Find tangent and normal
  1024. tangent = point(x1[idx_m], y1[idx_m]) - point(x1[idx_p], y1[idx_p])
  1025. tangent = unit_vec(tangent)
  1026. normal = point(tangent[1], -tangent[0])
  1027. # Displace length "thickness" in "positive" normal direction
  1028. displaced_pt = point(x1[idx], y1[idx]) + thickness*normal
  1029. x2[idx], y2[idx] = displaced_pt
  1030. for i in range(1, len(x1)-1):
  1031. displace(i-1, i+1, i) # centered difference for normal comp.
  1032. # One-sided differences at the end points
  1033. i = 0
  1034. displace(i, i+1, i)
  1035. i = len(x1)-1
  1036. displace(i-1, i, i)
  1037. # Combine x1,y1 with x2,y2 reversed
  1038. from numpy import concatenate
  1039. x = concatenate((x1, x2[-1::-1]))
  1040. y = concatenate((y1, y2[-1::-1]))
  1041. wall = Curve(x, y)
  1042. wall.set_filled_curves(color='white', pattern=pattern)
  1043. x = [x1[-1]] + x2[-1::-1].tolist() + [x1[0]]
  1044. y = [y1[-1]] + y2[-1::-1].tolist() + [y1[0]]
  1045. self.shapes['wall'] = wall
  1046. def geometric_features(self):
  1047. d = {'start': point(self.x1[0], self.y1[0]),
  1048. 'end': point(self.x1[-1], self.y1[-1])}
  1049. return d
  1050. class VelocityProfile(Shape):
  1051. def __init__(self, start, height, profile, num_arrows, scaling=1):
  1052. # vx, vy = profile(y)
  1053. shapes = {}
  1054. # Draw left line
  1055. shapes['start line'] = Line(start, (start[0], start[1]+height))
  1056. # Draw velocity arrows
  1057. dy = float(height)/(num_arrows-1)
  1058. x = start[0]
  1059. y = start[1]
  1060. r = profile(y) # Test on return type
  1061. if not isinstance(r, (list,tuple,ndarray)) and len(r) != 2:
  1062. raise TypeError('VelocityProfile constructor: profile(y) function must return velocity vector (vx,vy), not %s' % type(r))
  1063. for i in range(num_arrows):
  1064. y = start[1] + i*dy
  1065. vx, vy = profile(y)
  1066. if abs(vx) < 1E-8:
  1067. continue
  1068. vx *= scaling
  1069. vy *= scaling
  1070. arr = Arrow1((x,y), (x+vx, y+vy), '->')
  1071. shapes['arrow%d' % i] = arr
  1072. # Draw smooth profile
  1073. xs = []
  1074. ys = []
  1075. n = 100
  1076. dy = float(height)/n
  1077. for i in range(n+2):
  1078. y = start[1] + i*dy
  1079. vx, vy = profile(y)
  1080. vx *= scaling
  1081. vy *= scaling
  1082. xs.append(x+vx)
  1083. ys.append(y+vy)
  1084. shapes['smooth curve'] = Curve(xs, ys)
  1085. self.shapes = shapes
  1086. class Arrow1(Shape):
  1087. """Draw an arrow as Line with arrow."""
  1088. def __init__(self, start, end, style='->'):
  1089. arrow = Line(start, end)
  1090. arrow.set_arrow(style)
  1091. # Note:
  1092. self.shapes = {'arrow': arrow}
  1093. def geometric_features(self):
  1094. return self.shapes['arrow'].geometric_features()
  1095. class Arrow3(Shape):
  1096. """
  1097. Build a vertical line and arrow head from Line objects.
  1098. Then rotate `rotation_angle`.
  1099. """
  1100. def __init__(self, start, length, rotation_angle=0):
  1101. self.bottom = start
  1102. self.length = length
  1103. self.angle = rotation_angle
  1104. top = (self.bottom[0], self.bottom[1] + self.length)
  1105. main = Line(self.bottom, top)
  1106. #head_length = self.length/8.0
  1107. head_length = drawing_tool.xrange/50.
  1108. head_degrees = radians(30)
  1109. head_left_pt = (top[0] - head_length*sin(head_degrees),
  1110. top[1] - head_length*cos(head_degrees))
  1111. head_right_pt = (top[0] + head_length*sin(head_degrees),
  1112. top[1] - head_length*cos(head_degrees))
  1113. head_left = Line(head_left_pt, top)
  1114. head_right = Line(head_right_pt, top)
  1115. head_left.set_linestyle('solid')
  1116. head_right.set_linestyle('solid')
  1117. self.shapes = {'line': main, 'head left': head_left,
  1118. 'head right': head_right}
  1119. # rotate goes through self.shapes so self.shapes
  1120. # must be initialized first
  1121. self.rotate(rotation_angle, start)
  1122. def geometric_features(self):
  1123. return self.shapes['line'].geometric_features()
  1124. class Text(Point):
  1125. """
  1126. Place `text` at the (x,y) point `position`, with the given
  1127. fontsize (0 indicates that the default fontsize set in drawing_tool
  1128. is to be used). The text is centered around `position` if `alignment` is
  1129. 'center'; if 'left', the text starts at `position`, and if
  1130. 'right', the right and of the text is located at `position`.
  1131. """
  1132. def __init__(self, text, position, alignment='center', fontsize=0):
  1133. is_sequence(position)
  1134. is_sequence(position, length=2, can_be_None=True)
  1135. self.text = text
  1136. self.position = position
  1137. self.alignment = alignment
  1138. self.fontsize = fontsize
  1139. Point.__init__(self, position[0], position[1])
  1140. #no need for self.shapes here
  1141. def draw(self):
  1142. drawing_tool.text(self.text, (self.x, self.y),
  1143. self.alignment, self.fontsize)
  1144. def __str__(self):
  1145. return 'text "%s" at (%g,%g)' % (self.text, self.x, self.y)
  1146. def __repr__(self):
  1147. return repr(str(self))
  1148. class Text_wArrow(Text):
  1149. """
  1150. As class Text, but an arrow is drawn from the mid part of the text
  1151. to some point `arrow_tip`.
  1152. """
  1153. def __init__(self, text, position, arrow_tip,
  1154. alignment='center', fontsize=0):
  1155. is_sequence(arrow_tip, length=2, can_be_None=True)
  1156. is_sequence(position)
  1157. self.arrow_tip = arrow_tip
  1158. Text.__init__(self, text, position, alignment, fontsize)
  1159. def draw(self):
  1160. drawing_tool.text(self.text, self.position,
  1161. self.alignment, self.fontsize,
  1162. self.arrow_tip)
  1163. def __str__(self):
  1164. return 'annotation "%s" at (%g,%g) with arrow to (%g,%g)' % \
  1165. (self.text, self.x, self.y,
  1166. self.arrow_tip[0], self.arrow_tip[1])
  1167. def __repr__(self):
  1168. return repr(str(self))
  1169. class Axis(Shape):
  1170. def __init__(self, start, length, label,
  1171. rotation_angle=0, fontsize=0,
  1172. label_spacing=1./45, label_alignment='left'):
  1173. """
  1174. Draw axis from start with `length` to the right
  1175. (x axis). Place label at the end of the arrow tip.
  1176. Then return `rotation_angle` (in degrees).
  1177. The `label_spacing` denotes the space between the label
  1178. and the arrow tip as a fraction of the length of the plot
  1179. in x direction. A tuple can be given to adjust the position
  1180. in both the x and y directions (with one parameter, the
  1181. x position is adjusted).
  1182. With `label_alignment` one can place
  1183. the axis label text such that the arrow tip is to the 'left',
  1184. 'right', or 'center' with respect to the text field.
  1185. The `label_spacing` and `label_alignment`parameters can
  1186. be used to fine-tune the location of the label.
  1187. """
  1188. # Arrow is vertical arrow, make it horizontal
  1189. arrow = Arrow3(start, length, rotation_angle=-90)
  1190. arrow.rotate(rotation_angle, start)
  1191. if isinstance(label_spacing, (list,tuple)) and len(label_spacing) == 2:
  1192. x_spacing = drawing_tool.xrange*label_spacing[0]
  1193. y_spacing = drawing_tool.yrange*label_spacing[1]
  1194. elif isinstance(label_spacing, (int,float)):
  1195. # just x spacing
  1196. x_spacing = drawing_tool.xrange*label_spacing
  1197. y_spacing = 0
  1198. # should increase spacing for downward pointing axis
  1199. label_pos = [start[0] + length + x_spacing, start[1] + y_spacing]
  1200. label = Text(label, position=label_pos, fontsize=fontsize)
  1201. label.rotate(rotation_angle, start)
  1202. self.shapes = {'arrow': arrow, 'label': label}
  1203. def geometric_features(self):
  1204. return self.shapes['arrow'].geometric_features()
  1205. # Maybe Axis3 with label below/above?
  1206. class Force(Arrow1):
  1207. """
  1208. Indication of a force by an arrow and a text (symbol). Draw an
  1209. arrow, starting at `start` and with the tip at `end`. The symbol
  1210. is placed at `text_pos`, which can be 'start', 'end' or the
  1211. coordinates of a point. If 'end' or 'start', the text is placed at
  1212. a distance `text_spacing` times the width of the total plotting
  1213. area away from the specified point.
  1214. """
  1215. def __init__(self, start, end, text, text_spacing=1./60,
  1216. fontsize=0, text_pos='start', text_alignment='center'):
  1217. Arrow1.__init__(self, start, end, style='->')
  1218. if isinstance(text_spacing, (tuple,list)):
  1219. if len(text_spacing) == 2:
  1220. spacing = point(drawing_tool.xrange*text_spacing[0],
  1221. drawing_tool.xrange*text_spacing[1])
  1222. else:
  1223. spacing = drawing_tool.xrange*text_spacing[0]
  1224. else:
  1225. # just a number, this is x spacing
  1226. spacing = drawing_tool.xrange*text_spacing
  1227. start, end = arr2D(start), arr2D(end)
  1228. # Two cases: label at bottom of line or top, need more
  1229. # spacing if bottom
  1230. downward = (end-start)[1] < 0
  1231. upward = not downward # for easy code reading
  1232. if isinstance(text_pos, str):
  1233. if text_pos == 'start':
  1234. spacing_dir = unit_vec(start - end)
  1235. if upward:
  1236. spacing *= 1.7
  1237. if isinstance(spacing, (int, float)):
  1238. text_pos = start + spacing*spacing_dir
  1239. else:
  1240. text_pos = start + spacing
  1241. elif text_pos == 'end':
  1242. spacing_dir = unit_vec(end - start)
  1243. if downward:
  1244. spacing *= 1.7
  1245. if isinstance(spacing, (int, float)):
  1246. text_pos = end + spacing*spacing_dir
  1247. else:
  1248. text_pos = end + spacing
  1249. self.shapes['text'] = Text(text, text_pos, fontsize=fontsize,
  1250. alignment=text_alignment)
  1251. def geometric_features(self):
  1252. d = Arrow1.geometric_features(self)
  1253. d['symbol_location'] = self.shapes['text'].position
  1254. return d
  1255. class Axis2(Force):
  1256. def __init__(self, start, length, label,
  1257. rotation_angle=0, fontsize=0,
  1258. label_spacing=1./45, label_alignment='left'):
  1259. direction = point(cos(radians(rotation_angle)),
  1260. sin(radians(rotation_angle)))
  1261. Force.__init__(start=start, end=length*direction, text=label,
  1262. text_spacing=label_spacing,
  1263. fontsize=fontsize, text_pos='end',
  1264. text_alignment=label_alignment)
  1265. # Substitute text by label for axis
  1266. self.shapes['label'] = self.shapes['text']
  1267. del self.shapes['text']
  1268. # geometric features from Force is ok
  1269. class Gravity(Axis):
  1270. """Downward-pointing gravity arrow with the symbol g."""
  1271. def __init__(self, start, length, fontsize=0):
  1272. Axis.__init__(self, start, length, '$g$', below=False,
  1273. rotation_angle=-90, label_spacing=1./30,
  1274. fontsize=fontsize)
  1275. self.shapes['arrow'].set_linecolor('black')
  1276. class Gravity(Force):
  1277. """Downward-pointing gravity arrow with the symbol g."""
  1278. def __init__(self, start, length, text='$g$', fontsize=0):
  1279. Force.__init__(self, start, (start[0], start[1]-length),
  1280. text, text_spacing=1./60,
  1281. fontsize=0, text_pos='end')
  1282. self.shapes['arrow'].set_linecolor('black')
  1283. class Distance_wText(Shape):
  1284. """
  1285. Arrow <-> with text (usually a symbol) at the midpoint, used for
  1286. identifying a some distance in a figure. The text is placed
  1287. slightly to the right of vertical-like arrows, with text displaced
  1288. `text_spacing` times to total distance in x direction of the plot
  1289. area. The text is by default aligned 'left' in this case. For
  1290. horizontal-like arrows, the text is placed the same distance
  1291. above, but aligned 'center' by default (when `alignment` is None).
  1292. """
  1293. def __init__(self, start, end, text, fontsize=0, text_spacing=1/60.,
  1294. alignment=None, text_pos='mid'):
  1295. start = arr2D(start)
  1296. end = arr2D(end)
  1297. # Decide first if we have a vertical or horizontal arrow
  1298. vertical = abs(end[0]-start[0]) < 2*abs(end[1]-start[1])
  1299. if vertical:
  1300. # Assume end above start
  1301. if end[1] < start[1]:
  1302. start, end = end, start
  1303. if alignment is None:
  1304. alignment = 'left'
  1305. else: # horizontal arrow
  1306. # Assume start to the right of end
  1307. if start[0] < end[0]:
  1308. start, end = end, start
  1309. if alignment is None:
  1310. alignment = 'center'
  1311. tangent = end - start
  1312. # Tangeng goes always to the left and upward
  1313. normal = unit_vec([tangent[1], -tangent[0]])
  1314. mid = 0.5*(start + end) # midpoint of start-end line
  1315. if text_pos == 'mid':
  1316. text_pos = mid + normal*drawing_tool.xrange*text_spacing
  1317. text = Text(text, text_pos, fontsize=fontsize,
  1318. alignment=alignment)
  1319. else:
  1320. is_sequence(text_pos, length=2)
  1321. text = Text_wArrow(text, text_pos, mid, alignment='left',
  1322. fontsize=fontsize)
  1323. arrow = Arrow1(start, end, style='<->')
  1324. arrow.set_linecolor('black')
  1325. arrow.set_linewidth(1)
  1326. self.shapes = {'arrow': arrow, 'text': text}
  1327. def geometric_features(self):
  1328. d = self.shapes['arrow'].geometric_features()
  1329. d['text_position'] = self.shapes['text'].position
  1330. return d
  1331. class Arc_wText(Shape):
  1332. def __init__(self, text, center, radius,
  1333. start_angle, arc_angle, fontsize=0,
  1334. resolution=180, text_spacing=1/60.):
  1335. arc = Arc(center, radius, start_angle, arc_angle,
  1336. resolution)
  1337. mid = arr2D(arc(arc_angle/2.))
  1338. normal = unit_vec(mid - arr2D(center))
  1339. text_pos = mid + normal*drawing_tool.xrange*text_spacing
  1340. self.shapes = {'arc': arc,
  1341. 'text': Text(text, text_pos, fontsize=fontsize)}
  1342. class Composition(Shape):
  1343. def __init__(self, shapes):
  1344. """shapes: list or dict of Shape objects."""
  1345. if isinstance(shapes, (tuple,list)):
  1346. # Convert to dict using the type of the list element as key
  1347. # (add a counter to make the keys unique)
  1348. shapes = {s.__class__.__name__ + '_' + str(i): s
  1349. for i, s in enumerate(shapes)}
  1350. self.shapes = shapes
  1351. # can make help methods: Line.midpoint, Line.normal(pt, dir='left') -> (x,y)
  1352. # list annotations in each class? contains extra annotations for explaining
  1353. # important parameters to the constructor, e.g., Line.annotations holds
  1354. # start and end as Text objects. Shape.demo calls shape.draw and
  1355. # for annotation in self.demo: annotation.draw() YES!
  1356. # Can make overall demo of classes by making objects and calling demo
  1357. # Could include demo fig in each constructor
  1358. class SimplySupportedBeam(Shape):
  1359. def __init__(self, pos, size):
  1360. pos = arr2D(pos)
  1361. P0 = (pos[0] - size/2., pos[1]-size)
  1362. P1 = (pos[0] + size/2., pos[1]-size)
  1363. triangle = Triangle(P0, P1, pos)
  1364. gap = size/5.
  1365. h = size/4. # height of rectangle
  1366. P2 = (P0[0], P0[1]-gap-h)
  1367. rectangle = Rectangle(P2, size, h).set_filled_curves(pattern='/')
  1368. self.shapes = {'triangle': triangle, 'rectangle': rectangle}
  1369. self.dimensions = {'pos': Text('pos', pos),
  1370. 'size': Distance_wText((P2[0], P2[1]-size),
  1371. (P2[0]+size, P2[1]-size),
  1372. 'size')}
  1373. def geometric_features(self):
  1374. t = self.shapes['triangle']
  1375. r = self.shapes['rectangle']
  1376. d = {'pos': t.geometric_features()['p2'],
  1377. 'mid_support': r.geometric_features()['lower_mid']}
  1378. return d
  1379. class ConstantBeamLoad(Shape):
  1380. """
  1381. Downward-pointing arrows indicating a vertical load.
  1382. The arrows are of equal length and filling a rectangle
  1383. specified as in the :class:`Rectangle` class.
  1384. Recorded geometric features:
  1385. ==================== =============================================
  1386. Attribute Description
  1387. ==================== =============================================
  1388. mid_top Middle point at the top of the row of
  1389. arrows (often used for positioning a text).
  1390. ==================== =============================================
  1391. """
  1392. def __init__(self, lower_left_corner, width, height, num_arrows=10):
  1393. box = Rectangle(lower_left_corner, width, height)
  1394. self.shapes = {'box': box}
  1395. dx = float(width)/(num_arrows-1)
  1396. y_top = lower_left_corner[1] + height
  1397. y_tip = lower_left_corner[1]
  1398. for i in range(num_arrows):
  1399. x = lower_left_corner[0] + i*dx
  1400. self.shapes['arrow%d' % i] = Arrow1((x, y_top), (x, y_tip))
  1401. def geometric_features(self):
  1402. return {'mid_top': self.shapes['box'].geometric_features()['upper_mid']}
  1403. class Moment(Arc_wText):
  1404. def __init__(self, text, center, radius,
  1405. left=True, counter_clockwise=True,
  1406. fontsize=0, text_spacing=1/60.):
  1407. style = '->' if counter_clockwise else '<-'
  1408. start_angle = 90 if left else -90
  1409. Arc_wText.__init__(self, text, center, radius,
  1410. start_angle=start_angle,
  1411. arc_angle=180, fontsize=fontsize,
  1412. text_spacing=text_spacing,
  1413. resolution=180)
  1414. self.shapes['arc']['arc'].set_arrow(style) # Curve object
  1415. class Wheel(Shape):
  1416. def __init__(self, center, radius, inner_radius=None, nlines=10):
  1417. if inner_radius is None:
  1418. inner_radius = radius/5.0
  1419. outer = Circle(center, radius)
  1420. inner = Circle(center, inner_radius)
  1421. lines = []
  1422. # Draw nlines+1 since the first and last coincide
  1423. # (then nlines lines will be visible)
  1424. t = linspace(0, 2*pi, self.nlines+1)
  1425. Ri = inner_radius; Ro = radius
  1426. x0 = center[0]; y0 = center[1]
  1427. xinner = x0 + Ri*cos(t)
  1428. yinner = y0 + Ri*sin(t)
  1429. xouter = x0 + Ro*cos(t)
  1430. youter = y0 + Ro*sin(t)
  1431. lines = [Line((xi,yi),(xo,yo)) for xi, yi, xo, yo in \
  1432. zip(xinner, yinner, xouter, youter)]
  1433. self.shapes = {'inner': inner, 'outer': outer,
  1434. 'spokes': Composition(
  1435. {'spoke%d' % i: lines[i]
  1436. for i in range(len(lines))})}
  1437. class SineWave(Shape):
  1438. def __init__(self, xstart, xstop,
  1439. wavelength, amplitude, mean_level):
  1440. self.xstart = xstart
  1441. self.xstop = xstop
  1442. self.wavelength = wavelength
  1443. self.amplitude = amplitude
  1444. self.mean_level = mean_level
  1445. npoints = (self.xstop - self.xstart)/(self.wavelength/61.0)
  1446. x = linspace(self.xstart, self.xstop, npoints)
  1447. k = 2*pi/self.wavelength # frequency
  1448. y = self.mean_level + self.amplitude*sin(k*x)
  1449. self.shapes = {'waves': Curve(x,y)}
  1450. class Spring(Shape):
  1451. """
  1452. Specify a *vertical* spring, starting at `start` and with `length`
  1453. as total vertical length. In the middle of the spring there are
  1454. `num_windings` circular windings to illustrate the spring. If
  1455. `teeth` is true, the spring windings look like saw teeth,
  1456. otherwise the windings are smooth circles. The parameters `width`
  1457. (total width of spring) and `bar_length` (length of first and last
  1458. bar are given sensible default values if they are not specified
  1459. (these parameters can later be extracted as attributes, see table
  1460. below).
  1461. """
  1462. spring_fraction = 1./2 # fraction of total length occupied by spring
  1463. def __init__(self, start, length, width=None, bar_length=None,
  1464. num_windings=11, teeth=False):
  1465. B = start
  1466. n = num_windings - 1 # n counts teeth intervals
  1467. if n <= 6:
  1468. n = 7
  1469. # n must be odd:
  1470. if n % 2 == 0:
  1471. n = n+1
  1472. L = length
  1473. if width is None:
  1474. w = L/10.
  1475. else:
  1476. w = width/2.0
  1477. s = bar_length
  1478. # [0, x, L-x, L], f = (L-2*x)/L
  1479. # x = L*(1-f)/2.
  1480. # B: start point
  1481. # w: half-width
  1482. # L: total length
  1483. # s: length of first bar
  1484. # P0: start of dashpot (B[0]+s)
  1485. # P1: end of dashpot
  1486. # P2: end point
  1487. shapes = {}
  1488. if s is None:
  1489. f = Spring.spring_fraction
  1490. s = L*(1-f)/2. # start of spring
  1491. self.bar_length = s # record
  1492. self.width = 2*w
  1493. P0 = (B[0], B[1] + s)
  1494. P1 = (B[0], B[1] + L-s)
  1495. P2 = (B[0], B[1] + L)
  1496. if s >= L:
  1497. raise ValueError('length of first bar: %g is larger than total length: %g' % (s, L))
  1498. shapes['bar1'] = Line(B, P0)
  1499. spring_length = L - 2*s
  1500. t = spring_length/n # height increment per winding
  1501. if teeth:
  1502. resolution = 4
  1503. else:
  1504. resolution = 90
  1505. q = linspace(0, n, n*resolution + 1)
  1506. x = P0[0] + w*sin(2*pi*q)
  1507. y = P0[1] + q*t
  1508. shapes['spiral'] = Curve(x, y)
  1509. shapes['bar2'] = Line(P1,P2)
  1510. self.shapes = shapes
  1511. # Dimensions
  1512. start = Text_wArrow('start', (B[0]-1.5*w,B[1]-1.5*w), B)
  1513. width = Distance_wText((B[0]-w, B[1]-3.5*w), (B[0]+w, B[1]-3.5*w),
  1514. 'width')
  1515. length = Distance_wText((B[0]+3*w, B[1]), (B[0]+3*w, B[1]+L),
  1516. 'length')
  1517. num_windings = Text_wArrow('num_windings',
  1518. (B[0]+2*w,P2[1]+w),
  1519. (B[0]+1.2*w, B[1]+L/2.))
  1520. blength1 = Distance_wText((B[0]-2*w, B[1]), (B[0]-2*w, P0[1]),
  1521. 'bar_length',
  1522. text_pos=(P0[0]-7*w, P0[1]+w))
  1523. blength2 = Distance_wText((P1[0]-2*w, P1[1]), (P2[0]-2*w, P2[1]),
  1524. 'bar_length',
  1525. text_pos=(P2[0]-7*w, P2[1]+w))
  1526. dims = {'start': start, 'width': width, 'length': length,
  1527. 'num_windings': num_windings, 'bar_length1': blength1,
  1528. 'bar_length2': blength2}
  1529. self.dimensions = dims
  1530. def geometric_features(self):
  1531. """
  1532. Recorded geometric features:
  1533. ==================== =============================================
  1534. Attribute Description
  1535. ==================== =============================================
  1536. start Start point of spring.
  1537. end End point of spring.
  1538. width Total width of spring.
  1539. bar_length Length of first (and last) bar part.
  1540. ==================== =============================================
  1541. """
  1542. b1 = self.shapes['bar1']
  1543. d = {'start': b1.geometric_features()['start'],
  1544. 'end': self.shapes['bar2'].geometric_features()['end'],
  1545. 'bar_length': self.bar_length,
  1546. 'width': self.width}
  1547. return d
  1548. class Dashpot(Shape):
  1549. """
  1550. Specify a vertical dashpot of height `total_length` and `start` as
  1551. bottom/starting point. The first bar part has length `bar_length`.
  1552. Then comes the dashpot as a rectangular construction of total
  1553. width `width` and height `dashpot_length`. The position of the
  1554. piston inside the rectangular dashpot area is given by
  1555. `piston_pos`, which is the distance between the first bar (given
  1556. by `bar_length`) to the piston.
  1557. If some of `dashpot_length`, `bar_length`, `width` or `piston_pos`
  1558. are not given, suitable default values are calculated. Their
  1559. values can be extracted as keys in the dict returned from
  1560. ``geometric_features``.
  1561. """
  1562. dashpot_fraction = 1./2 # fraction of total_length
  1563. piston_gap_fraction = 1./6 # fraction of width
  1564. piston_thickness_fraction = 1./8 # fraction of dashplot_length
  1565. def __init__(self, start, total_length, bar_length=None,
  1566. width=None, dashpot_length=None, piston_pos=None):
  1567. B = start
  1568. L = total_length
  1569. if width is None:
  1570. w = L/10. # total width 1/5 of length
  1571. else:
  1572. w = width/2.0
  1573. s = bar_length
  1574. # [0, x, L-x, L], f = (L-2*x)/L
  1575. # x = L*(1-f)/2.
  1576. # B: start point
  1577. # w: half-width
  1578. # L: total length
  1579. # s: length of first bar
  1580. # P0: start of dashpot (B[0]+s)
  1581. # P1: end of dashpot
  1582. # P2: end point
  1583. shapes = {}
  1584. # dashpot is P0-P1 in y and width 2*w
  1585. if dashpot_length is None:
  1586. if s is None:
  1587. f = Dashpot.dashpot_fraction
  1588. s = L*(1-f)/2. # default
  1589. P1 = (B[0], B[1]+L-s)
  1590. dashpot_length = f*L
  1591. else:
  1592. if s is None:
  1593. f = 1./2 # the bar lengths are taken as f*dashpot_length
  1594. s = f*dashpot_length # default
  1595. P1 = (B[0], B[1]+s+dashpot_length)
  1596. P0 = (B[0], B[1]+s)
  1597. P2 = (B[0], B[1]+L)
  1598. if P2[1] > P1[1] > P0[1]:
  1599. pass # ok
  1600. else:
  1601. raise ValueError('Dashpot has inconsistent dimensions! start: %g, dashpot begin: %g, dashpot end: %g, very end: %g' % (B[1], P0[1], P1[1], P2[1]))
  1602. shapes['line start'] = Line(B, P0)
  1603. shapes['pot'] = Curve([P1[0]-w, P0[0]-w, P0[0]+w, P1[0]+w],
  1604. [P1[1], P0[1], P0[1], P1[1]])
  1605. piston_thickness = dashpot_length*Dashpot.piston_thickness_fraction
  1606. if piston_pos is None:
  1607. piston_pos = 1/3.*dashpot_length
  1608. if piston_pos < 0:
  1609. piston_pos = 0
  1610. elif piston_pos > dashpot_length:
  1611. piston_pos = dashpot_length - piston_thickness
  1612. abs_piston_pos = P0[1] + piston_pos
  1613. gap = w*Dashpot.piston_gap_fraction
  1614. shapes['piston'] = Composition(
  1615. {'line': Line(P2, (B[0], abs_piston_pos + piston_thickness)),
  1616. 'rectangle': Rectangle((B[0] - w+gap, abs_piston_pos),
  1617. 2*w-2*gap, piston_thickness),
  1618. })
  1619. shapes['piston']['rectangle'].set_filled_curves(pattern='X')
  1620. self.shapes = shapes
  1621. self.bar_length = s
  1622. self.width = 2*w
  1623. self.piston_pos = piston_pos
  1624. self.dashpot_length = dashpot_length
  1625. # Dimensions
  1626. start = Text_wArrow('start', (B[0]-1.5*w,B[1]-1.5*w), B)
  1627. width = Distance_wText((B[0]-w, B[1]-3.5*w), (B[0]+w, B[1]-3.5*w),
  1628. 'width')
  1629. dplength = Distance_wText((B[0]+2*w, P0[1]), (B[0]+2*w, P1[1]),
  1630. 'dashpot_length', text_pos=(B[0]+w,B[1]-w))
  1631. blength = Distance_wText((B[0]-2*w, B[1]), (B[0]-2*w, P0[1]),
  1632. 'bar_length', text_pos=(B[0]-6*w,P0[1]-w))
  1633. ppos = Distance_wText((B[0]-2*w, P0[1]), (B[0]-2*w, P0[1]+piston_pos),
  1634. 'piston_pos', text_pos=(B[0]-6*w,P0[1]+piston_pos-w))
  1635. tlength = Distance_wText((B[0]+4*w, B[1]), (B[0]+4*w, B[1]+L),
  1636. 'total_length',
  1637. text_pos=(B[0]+4.5*w, B[1]+L-2*w))
  1638. line = Line((B[0]+w, abs_piston_pos), (B[0]+7*w, abs_piston_pos)).set_linestyle('dashed').set_linecolor('black').set_linewidth(1)
  1639. pp = Text('abs_piston_pos', (B[0]+7*w, abs_piston_pos), alignment='left')
  1640. dims = {'start': start, 'width': width, 'dashpot_length': dplength,
  1641. 'bar_length': blength, 'total_length': tlength,
  1642. 'piston_pos': ppos,}
  1643. #'abs_piston_pos': Composition({'line': line, 'text': pp})}
  1644. self.dimensions = dims
  1645. def geometric_features(self):
  1646. """
  1647. Recorded geometric features:
  1648. ==================== =============================================
  1649. Attribute Description
  1650. ==================== =============================================
  1651. start Start point of dashpot.
  1652. end End point of dashpot.
  1653. bar_length Length of first bar (from start to spring).
  1654. dashpot_length Length of dashpot middle part.
  1655. width Total width of dashpot.
  1656. piston_pos Position of piston in dashpot, relative to
  1657. start[1] + bar_length.
  1658. ==================== =============================================
  1659. """
  1660. d = {'start': self.shapes['line start'].geometric_features()['start'],
  1661. 'end': self.shapes['piston']['line'].geometric_features()['start'],
  1662. 'bar_length': self.bar_length,
  1663. 'piston_pos': self.piston_pos,
  1664. 'width': self.width,
  1665. 'dashpot_length': self.dashpot_length,
  1666. }
  1667. return d
  1668. class Wavy(Shape):
  1669. """
  1670. A wavy graph consisting of a user-given main curve y=f(x) with
  1671. additional sinusoidal waves of given (constant) amplitude,
  1672. but varying wavelength (a characteristic wavelength is specified).
  1673. """
  1674. def __init__(self, main_curve, interval, wavelength_of_perturbations,
  1675. amplitude_of_perturbations, smoothness):
  1676. """
  1677. ============================ ====================================
  1678. Name Description
  1679. ============================ ====================================
  1680. main_curve f(x) Python function
  1681. interval interval for main_curve
  1682. wavelength_of_perturbations dominant wavelength perturbed waves
  1683. amplitude_of_perturbations amplitude of perturbed waves
  1684. smoothness in [0, 1]: smooth=0, rough=1
  1685. ============================ ====================================
  1686. """
  1687. xmin, xmax = interval
  1688. L = wavelength_of_perturbations
  1689. k_0 = 2*pi/L # main frequency of waves
  1690. k_p = k_0*0.5
  1691. k_k = k_0/2*smoothness
  1692. A_0 = amplitude_of_perturbations
  1693. A_p = 0.3*A_0
  1694. A_k = k_0/2
  1695. x = linspace(xmin, xmax, 2001)
  1696. def w(x):
  1697. A = A_0 + A_p*sin(A_k*x)
  1698. k = k_0 + k_p*sin(k_k*x)
  1699. y = main_curve(x) + A*sin(k*x)
  1700. return y
  1701. self.shapes = {'wavy': Curve(x, w(x))}
  1702. # Use closure w to define __call__ - then we do not need
  1703. # to store all the parameters A_0, A_k, etc. as attributes
  1704. self.__call__ = w
  1705. class StochasticWavyCurve(object):
  1706. """
  1707. Precomputed stochastic wavy graphs.
  1708. There are three graphs with different look.
  1709. Curve 0:
  1710. ----------------------------------------------------------------------
  1711. |
  1712. |
  1713. *|
  1714. * |
  1715. * |
  1716. * |
  1717. * |
  1718. * |
  1719. * |
  1720. * |
  1721. * |
  1722. * |
  1723. * |
  1724. * |
  1725. |*
  1726. | *
  1727. | *
  1728. | *
  1729. | *
  1730. | *
  1731. | *
  1732. | *
  1733. | *
  1734. | *
  1735. | *
  1736. | *
  1737. | *
  1738. | *
  1739. | *
  1740. | *
  1741. | *
  1742. | *
  1743. | *
  1744. | *
  1745. | *
  1746. | *
  1747. | *
  1748. | *
  1749. | *
  1750. | *
  1751. | *
  1752. | *
  1753. | *
  1754. | *
  1755. | *
  1756. | *
  1757. | *
  1758. | *
  1759. | *
  1760. | *
  1761. | *
  1762. | *
  1763. | *
  1764. | *
  1765. | *
  1766. | *
  1767. | *
  1768. | *
  1769. | *
  1770. | *
  1771. | *
  1772. | *
  1773. | *
  1774. | *
  1775. | *
  1776. | *
  1777. | *
  1778. | *
  1779. | *
  1780. | *
  1781. | *
  1782. | *
  1783. |*
  1784. *|
  1785. * |
  1786. * |
  1787. * |
  1788. * |
  1789. * |
  1790. * |
  1791. * |
  1792. * |
  1793. * |
  1794. * |
  1795. * |
  1796. * |
  1797. * |
  1798. * |
  1799. * |
  1800. |*
  1801. | *
  1802. | *
  1803. | *
  1804. | *
  1805. | *
  1806. | *
  1807. | *
  1808. | *
  1809. | *
  1810. | *
  1811. | *
  1812. | *
  1813. | *
  1814. | *
  1815. | *
  1816. | *
  1817. | *
  1818. | *
  1819. | *
  1820. | *
  1821. | *
  1822. | *
  1823. | *
  1824. | *
  1825. | *
  1826. | *
  1827. | *
  1828. | *
  1829. | *
  1830. | *
  1831. | *
  1832. | *
  1833. | *
  1834. | *
  1835. | *
  1836. | *
  1837. | *
  1838. | *
  1839. | *
  1840. | *
  1841. | *
  1842. | *
  1843. | *
  1844. | *
  1845. |
  1846. * |
  1847. * |
  1848. * |
  1849. * |
  1850. * |
  1851. * |
  1852. * |
  1853. * |
  1854. * |
  1855. * |
  1856. * |
  1857. * |
  1858. * |
  1859. * |
  1860. * |
  1861. * |
  1862. * |
  1863. * |
  1864. * |
  1865. * |
  1866. * |
  1867. * |
  1868. * |
  1869. * |
  1870. * |
  1871. * |
  1872. * |
  1873. * |
  1874. * |
  1875. * |
  1876. Curve 2:
  1877. ----------------------------------------------------------------------
  1878. |
  1879. |
  1880. |
  1881. |*
  1882. |*
  1883. |*
  1884. |
  1885. |
  1886. *|
  1887. |*
  1888. | *
  1889. | *
  1890. | *
  1891. | *
  1892. | *
  1893. | *
  1894. | *
  1895. | *
  1896. | *
  1897. | *
  1898. | *
  1899. | *
  1900. | *
  1901. | *
  1902. | *
  1903. | *
  1904. | *
  1905. | *
  1906. | *
  1907. | *
  1908. | *
  1909. | *
  1910. | *
  1911. | *
  1912. | *
  1913. | *
  1914. | *
  1915. | *
  1916. | *
  1917. | *
  1918. | *
  1919. | *
  1920. | *
  1921. | *
  1922. | *
  1923. | *
  1924. | *
  1925. | *
  1926. | *
  1927. | *
  1928. | *
  1929. | *
  1930. | *
  1931. | *
  1932. | *
  1933. | *
  1934. | *
  1935. | *
  1936. | *
  1937. | *
  1938. | *
  1939. |
  1940. * |
  1941. * |
  1942. * |
  1943. * |
  1944. * |
  1945. * |
  1946. * |
  1947. * |
  1948. * |
  1949. * |
  1950. * |
  1951. * |
  1952. * |
  1953. * |
  1954. * |
  1955. * |
  1956. * |
  1957. * |
  1958. * |
  1959. * |
  1960. * |
  1961. * |
  1962. * |
  1963. * |
  1964. * |
  1965. * |
  1966. * |
  1967. * |
  1968. * |
  1969. |
  1970. | *
  1971. | *
  1972. | *
  1973. | *
  1974. | *
  1975. | *
  1976. | *
  1977. | *
  1978. | *
  1979. | *
  1980. | *
  1981. | *
  1982. | *
  1983. | *
  1984. | *
  1985. | *
  1986. | *
  1987. | *
  1988. | *
  1989. | *
  1990. | *
  1991. | *
  1992. |*
  1993. |*
  1994. |
  1995. |
  1996. |
  1997. |*
  1998. | *
  1999. | *
  2000. |*
  2001. |
  2002. *|
  2003. |*
  2004. | *
  2005. | *
  2006. | *
  2007. | *
  2008. | *
  2009. | *
  2010. | *
  2011. | *
  2012. | *
  2013. | *
  2014. | *
  2015. | *
  2016. | *
  2017. | *
  2018. | *
  2019. | *
  2020. | *
  2021. | *
  2022. | *
  2023. | *
  2024. | *
  2025. | *
  2026. | *
  2027. | *
  2028. | *
  2029. | *
  2030. | *
  2031. | *
  2032. | *
  2033. | *
  2034. | *
  2035. | *
  2036. | *
  2037. | *
  2038. | *
  2039. | *
  2040. | *
  2041. | *
  2042. | *
  2043. Curve 2:
  2044. ----------------------------------------------------------------------
  2045. |
  2046. |
  2047. |
  2048. |
  2049. |*
  2050. | *
  2051. | *
  2052. | *
  2053. | *
  2054. | *
  2055. | *
  2056. | *
  2057. | *
  2058. | *
  2059. | *
  2060. | *
  2061. | *
  2062. | *
  2063. | *
  2064. | *
  2065. | *
  2066. | *
  2067. | *
  2068. | *
  2069. | *
  2070. | *
  2071. | *
  2072. | *
  2073. |*
  2074. |
  2075. * |
  2076. * |
  2077. * |
  2078. * |
  2079. * |
  2080. * |
  2081. * |
  2082. * |
  2083. * |
  2084. * |
  2085. |*
  2086. | *
  2087. | *
  2088. | *
  2089. | *
  2090. | *
  2091. | *
  2092. | *
  2093. | *
  2094. | *
  2095. | *
  2096. | *
  2097. | *
  2098. | *
  2099. | *
  2100. | *
  2101. | *
  2102. | *
  2103. | *
  2104. | *
  2105. *|
  2106. * |
  2107. * |
  2108. * |
  2109. * |
  2110. * |
  2111. * |
  2112. * |
  2113. * |
  2114. * |
  2115. * |
  2116. * |
  2117. * |
  2118. * |
  2119. * |
  2120. * |
  2121. * |
  2122. |
  2123. | *
  2124. | *
  2125. | *
  2126. | *
  2127. | *
  2128. | *
  2129. | *
  2130. | *
  2131. | *
  2132. | *
  2133. | *
  2134. | *
  2135. | *
  2136. | *
  2137. | *
  2138. | *
  2139. | *
  2140. | *
  2141. | *
  2142. | *
  2143. | *
  2144. | *
  2145. |*
  2146. *|
  2147. * |
  2148. * |
  2149. * |
  2150. * |
  2151. * |
  2152. * |
  2153. * |
  2154. * |
  2155. * |
  2156. * |
  2157. * |
  2158. * |
  2159. * |
  2160. * |
  2161. * |
  2162. * |
  2163. * |
  2164. * |
  2165. * |
  2166. * |
  2167. * |
  2168. * |
  2169. * |
  2170. * |
  2171. * |
  2172. * |
  2173. * |
  2174. * |
  2175. * |
  2176. * |
  2177. * |
  2178. * |
  2179. * |
  2180. * |
  2181. * |
  2182. * |
  2183. * |
  2184. * |
  2185. * |
  2186. * |
  2187. * |
  2188. * |
  2189. * |
  2190. * |
  2191. * |
  2192. * |
  2193. * |
  2194. * |
  2195. * |
  2196. * |
  2197. * |
  2198. * |
  2199. * |
  2200. * |
  2201. * |
  2202. *|
  2203. |*
  2204. | *
  2205. | *
  2206. | *
  2207. | *
  2208. | *
  2209. | *
  2210. See also hplgit.github.io/pysketcher/doc/src/tut/fig-tut/StochasticWavyCurve.png (and .pdf)
  2211. """
  2212. # The curves were generated by the script generate_road_profiles.py and
  2213. # the code below were generated by plot_roads.py. Both scripts are
  2214. # found doc/src/src-bumpy in the repo git@github.com:hplgit/bumpy.git
  2215. def __init__(self, curve_no=0, percentage=100):
  2216. """
  2217. ============= ===================================================
  2218. Argument Explanation
  2219. ============= ===================================================
  2220. curve_no 0, 1, or 2: chooses one out of three shapes.
  2221. percentage The percentage of the defined curve to be used.
  2222. ============= ===================================================
  2223. """
  2224. self._define_curves()
  2225. self.curve_no = curve_no
  2226. m = int(len(self.x)/float(percentage)*100)
  2227. self.shapes = {'wavy': Curve(self.x[:m], self.y[curve_no][:m])}
  2228. def __call__(self, x):
  2229. raise NotImplementedError
  2230. def _define_curves(self):
  2231. self.x = array([0.0000, 0.0606, 0.1212, 0.1818, 0.2424, 0.3030, 0.3636, 0.4242, 0.4848, 0.5455, 0.6061, 0.6667, 0.7273, 0.7879, 0.8485, 0.9091, 0.9697, 1.0303, 1.0909, 1.1515, 1.2121, 1.2727, 1.3333, 1.3939, 1.4545, 1.5152, 1.5758, 1.6364, 1.6970, 1.7576, 1.8182, 1.8788, 1.9394, 2.0000, 2.0606, 2.1212, 2.1818, 2.2424, 2.3030, 2.3636, 2.4242, 2.4848, 2.5455, 2.6061, 2.6667, 2.7273, 2.7879, 2.8485, 2.9091, 2.9697, 3.0303, 3.0909, 3.1515, 3.2121, 3.2727, 3.3333, 3.3939, 3.4545, 3.5152, 3.5758, 3.6364, 3.6970, 3.7576, 3.8182, 3.8788, 3.9394, 4.0000, 4.0606, 4.1212, 4.1818, 4.2424, 4.3030, 4.3636, 4.4242, 4.4848, 4.5455, 4.6061, 4.6667, 4.7273, 4.7879, 4.8485, 4.9091, 4.9697, 5.0303, 5.0909, 5.1515, 5.2121, 5.2727, 5.3333, 5.3939, 5.4545, 5.5152, 5.5758, 5.6364, 5.6970, 5.7576, 5.8182, 5.8788, 5.9394, 6.0000, 6.0606, 6.1212, 6.1818, 6.2424, 6.3030, 6.3636, 6.4242, 6.4848, 6.5455, 6.6061, 6.6667, 6.7273, 6.7879, 6.8485, 6.9091, 6.9697, 7.0303, 7.0909, 7.1515, 7.2121, 7.2727, 7.3333, 7.3939, 7.4545, 7.5152, 7.5758, 7.6364, 7.6970, 7.7576, 7.8182, 7.8788, 7.9394, 8.0000, 8.0606, 8.1212, 8.1818, 8.2424, 8.3030, 8.3636, 8.4242, 8.4848, 8.5455, 8.6061, 8.6667, 8.7273, 8.7879, 8.8485, 8.9091, 8.9697, 9.0303, 9.0909, 9.1515, 9.2121, 9.2727, 9.3333, 9.3939, 9.4545, 9.5152, 9.5758, 9.6364, 9.6970, 9.7576, 9.8182, 9.8788, 9.9394, 10.0000, 10.0606, 10.1212, 10.1818, 10.2424, 10.3030, 10.3636, 10.4242, 10.4848, 10.5455, 10.6061, 10.6667, 10.7273, 10.7879, 10.8485, 10.9091, 10.9697, 11.0303, 11.0909, 11.1515, 11.2121, 11.2727, 11.3333, 11.3939, 11.4545, 11.5152, 11.5758, 11.6364, 11.6970, 11.7576, 11.8182, 11.8788, 11.9394, 12.0000, 12.0606, 12.1212, 12.1818, 12.2424, 12.3030, 12.3636, 12.4242, 12.4848, 12.5455, 12.6061, 12.6667, 12.7273, 12.7879, 12.8485, 12.9091, 12.9697, 13.0303, 13.0909, 13.1515, 13.2121, 13.2727, 13.3333, 13.3939, 13.4545, 13.5152, 13.5758, 13.6364, 13.6970, 13.7576, 13.8182, 13.8788, 13.9394, 14.0000, 14.0606, 14.1212, 14.1818, 14.2424, 14.3030, 14.3636, 14.4242, 14.4848, 14.5455, 14.6061, 14.6667, 14.7273, 14.7879, 14.8485, 14.9091, 14.9697, 15.0303, 15.0909, 15.1515, 15.2121, 15.2727, 15.3333, 15.3939, 15.4545, 15.5152, 15.5758, 15.6364, 15.6970, 15.7576, 15.8182, 15.8788, 15.9394, 16.0000, 16.0606, 16.1212, 16.1818, 16.2424, 16.3030, 16.3636, 16.4242, 16.4848, 16.5455, 16.6061, 16.6667, 16.7273, 16.7879, 16.8485, 16.9091, 16.9697, 17.0303, 17.0909, 17.1515, 17.2121, 17.2727, 17.3333, 17.3939, 17.4545, 17.5152, 17.5758, 17.6364, 17.6970, 17.7576, 17.8182, 17.8788, 17.9394, 18.0000, 18.0606, 18.1212, 18.1818, 18.2424, 18.3030, 18.3636, 18.4242, 18.4848, 18.5455, 18.6061, 18.6667, 18.7273, 18.7879, 18.8485, 18.9091, 18.9697, 19.0303, 19.0909, 19.1515, 19.2121, 19.2727, 19.3333, 19.3939, 19.4545, 19.5152, 19.5758, 19.6364, 19.6970, 19.7576, 19.8182, 19.8788, 19.9394, 20.0000, 20.0606, 20.1212, 20.1818, 20.2424, 20.3030, 20.3636, 20.4242, 20.4848, 20.5455, 20.6061, 20.6667, 20.7273, 20.7879, 20.8485, 20.9091, 20.9697, 21.0303, 21.0909, 21.1515, 21.2121, 21.2727, 21.3333, 21.3939, 21.4545, 21.5152, 21.5758, 21.6364, 21.6970, 21.7576, 21.8182, 21.8788, 21.9394, 22.0000, 22.0606, 22.1212, 22.1818, 22.2424, 22.3030, 22.3636, 22.4242, 22.4848, 22.5455, 22.6061, 22.6667, 22.7273, 22.7879, 22.8485, 22.9091, 22.9697, 23.0303, 23.0909, 23.1515, 23.2121, 23.2727, 23.3333, 23.3939, 23.4545, 23.5152, 23.5758, 23.6364, 23.6970, 23.7576, 23.8182, 23.8788, 23.9394, 24.0000, 24.0606, 24.1212, 24.1818, 24.2424, 24.3030, 24.3636, 24.4242, 24.4848, 24.5455, 24.6061, 24.6667, 24.7273, 24.7879, 24.8485, 24.9091, 24.9697, 25.0303, 25.0909, 25.1515, 25.2121, 25.2727, 25.3333, 25.3939, 25.4545, 25.5152, 25.5758, 25.6364, 25.6970, 25.7576, 25.8182, 25.8788, 25.9394, 26.0000, 26.0606, 26.1212, 26.1818, 26.2424, 26.3030, 26.3636, 26.4242, 26.4848, 26.5455, 26.6061, 26.6667, 26.7273, 26.7879, 26.8485, 26.9091, 26.9697, 27.0303, 27.0909, 27.1515, 27.2121, 27.2727, 27.3333, 27.3939, 27.4545, 27.5152, 27.5758, 27.6364, 27.6970, 27.7576, 27.8182, 27.8788, 27.9394, 28.0000, 28.0606, 28.1212, 28.1818, 28.2424, 28.3030, 28.3636, 28.4242, 28.4848, 28.5455, 28.6061, 28.6667, 28.7273, 28.7879, 28.8485, 28.9091, 28.9697, 29.0303, 29.0909, 29.1515, 29.2121, 29.2727, 29.3333, 29.3939, 29.4545, 29.5152, 29.5758, 29.6364, 29.6970, 29.7576, 29.8182, 29.8788, 29.9394, 30.0000, 30.0606, 30.1212, 30.1818, 30.2424, 30.3030, 30.3636, 30.4242, 30.4848, 30.5455, 30.6061, 30.6667, 30.7273, 30.7879, 30.8485, 30.9091, 30.9697, 31.0303, 31.0909, 31.1515, 31.2121, 31.2727, 31.3333, 31.3939, 31.4545, 31.5152, 31.5758, 31.6364, 31.6970, 31.7576, 31.8182, 31.8788, 31.9394, 32.0000, 32.0606, 32.1212, 32.1818, 32.2424, 32.3030, 32.3636, 32.4242, 32.4848, 32.5455, 32.6061, 32.6667, 32.7273, 32.7879, 32.8485, 32.9091, 32.9697, 33.0303, 33.0909, 33.1515, 33.2121, 33.2727, 33.3333, 33.3939, 33.4545, 33.5152, 33.5758, 33.6364, 33.6970, 33.7576, 33.8182, 33.8788, 33.9394, 34.0000, 34.0606, 34.1212, 34.1818, 34.2424, 34.3030, 34.3636, 34.4242, 34.4848, 34.5455, 34.6061, 34.6667, 34.7273, 34.7879, 34.8485, 34.9091, 34.9697, 35.0303, 35.0909, 35.1515, 35.2121, 35.2727, 35.3333, 35.3939, 35.4545, 35.5152, 35.5758, 35.6364, 35.6970, 35.7576, 35.8182, 35.8788, 35.9394, 36.0000, 36.0606, 36.1212, 36.1818, 36.2424, 36.3030, 36.3636, 36.4242, 36.4848, 36.5455, 36.6061, 36.6667, 36.7273, 36.7879, 36.8485, 36.9091, 36.9697, 37.0303, 37.0909, 37.1515, 37.2121, 37.2727, 37.3333, 37.3939, 37.4545, 37.5152, 37.5758, 37.6364, 37.6970, 37.7576, 37.8182, 37.8788, 37.9394, 38.0000, 38.0606, 38.1212, 38.1818, 38.2424, 38.3030, 38.3636, 38.4242, 38.4848, 38.5455, 38.6061, 38.6667, 38.7273, 38.7879, 38.8485, 38.9091, 38.9697, 39.0303, 39.0909, 39.1515, 39.2121, 39.2727, 39.3333, 39.3939, 39.4545, 39.5152, 39.5758, 39.6364, 39.6970, 39.7576, 39.8182, 39.8788, 39.9394, 40.0000, 40.0606, 40.1212, 40.1818, 40.2424, 40.3030, 40.3636, 40.4242, 40.4848, 40.5455, 40.6061, 40.6667, 40.7273, 40.7879, 40.8485, 40.9091, 40.9697, 41.0303, 41.0909, 41.1515, 41.2121, 41.2727, 41.3333, 41.3939, 41.4545, 41.5152, 41.5758, 41.6364, 41.6970, 41.7576, 41.8182, 41.8788, 41.9394, 42.0000, 42.0606, 42.1212, 42.1818, 42.2424, 42.3030, 42.3636, 42.4242, 42.4848, 42.5455, 42.6061, 42.6667, 42.7273, 42.7879, 42.8485, 42.9091, 42.9697, 43.0303, 43.0909, 43.1515, 43.2121, 43.2727, 43.3333, 43.3939, 43.4545, 43.5152, 43.5758, 43.6364, 43.6970, 43.7576, 43.8182, 43.8788, 43.9394, 44.0000, 44.0606, 44.1212, 44.1818, 44.2424, 44.3030, 44.3636, 44.4242, 44.4848, 44.5455, 44.6061, 44.6667, 44.7273, 44.7879, 44.8485, 44.9091, 44.9697, 45.0303, 45.0909, 45.1515, 45.2121, 45.2727, 45.3333, 45.3939, 45.4545, 45.5152, 45.5758, 45.6364, 45.6970, 45.7576, 45.8182, 45.8788, 45.9394, 46.0000, 46.0606, 46.1212, 46.1818, 46.2424, 46.3030, 46.3636, 46.4242, 46.4848, 46.5455, 46.6061, 46.6667, 46.7273, 46.7879, 46.8485, 46.9091, 46.9697, 47.0303, 47.0909, 47.1515, 47.2121, 47.2727, 47.3333, 47.3939, 47.4545, 47.5152, 47.5758, 47.6364, 47.6970, 47.7576, 47.8182, 47.8788, 47.9394, 48.0000, 48.0606, 48.1212, 48.1818, 48.2424, 48.3030, 48.3636, 48.4242, 48.4848, 48.5455, 48.6061, 48.6667, 48.7273, 48.7879, 48.8485, 48.9091, 48.9697, 49.0303, 49.0909, 49.1515, 49.2121, 49.2727, 49.3333, 49.3939, 49.4545, 49.5152, 49.5758, 49.6364, 49.6970, 49.7576, 49.8182, 49.8788, 49.9394, ])
  2232. self.y = [None]*3
  2233. self.y[0] = array([0.0000, 0.0005, 0.0006, 0.0004, -0.0004, -0.0007, -0.0022, -0.0027, -0.0036, -0.0042, -0.0050, -0.0049, -0.0060, -0.0072, -0.0085, -0.0092, -0.0104, -0.0116, -0.0133, -0.0148, -0.0160, -0.0177, -0.0186, -0.0191, -0.0192, -0.0187, -0.0187, -0.0187, -0.0192, -0.0198, -0.0201, -0.0208, -0.0216, -0.0227, -0.0242, -0.0260, -0.0277, -0.0299, -0.0319, -0.0328, -0.0333, -0.0338, -0.0347, -0.0360, -0.0363, -0.0365, -0.0370, -0.0373, -0.0364, -0.0355, -0.0343, -0.0329, -0.0317, -0.0312, -0.0309, -0.0306, -0.0301, -0.0290, -0.0275, -0.0259, -0.0238, -0.0222, -0.0200, -0.0176, -0.0154, -0.0130, -0.0108, -0.0081, -0.0046, -0.0001, 0.0035, 0.0061, 0.0083, 0.0105, 0.0130, 0.0156, 0.0170, 0.0181, 0.0196, 0.0212, 0.0231, 0.0247, 0.0262, 0.0277, 0.0293, 0.0309, 0.0325, 0.0336, 0.0348, 0.0360, 0.0378, 0.0401, 0.0423, 0.0443, 0.0457, 0.0473, 0.0488, 0.0500, 0.0511, 0.0518, 0.0528, 0.0534, 0.0547, 0.0561, 0.0577, 0.0585, 0.0594, 0.0606, 0.0611, 0.0614, 0.0617, 0.0612, 0.0607, 0.0608, 0.0603, 0.0599, 0.0588, 0.0577, 0.0557, 0.0543, 0.0532, 0.0520, 0.0505, 0.0496, 0.0499, 0.0490, 0.0489, 0.0496, 0.0504, 0.0504, 0.0509, 0.0512, 0.0512, 0.0504, 0.0499, 0.0498, 0.0493, 0.0491, 0.0483, 0.0478, 0.0474, 0.0468, 0.0462, 0.0460, 0.0462, 0.0467, 0.0472, 0.0476, 0.0483, 0.0491, 0.0502, 0.0510, 0.0504, 0.0503, 0.0514, 0.0527, 0.0538, 0.0547, 0.0554, 0.0561, 0.0561, 0.0558, 0.0548, 0.0540, 0.0531, 0.0524, 0.0516, 0.0513, 0.0511, 0.0520, 0.0519, 0.0513, 0.0512, 0.0525, 0.0535, 0.0545, 0.0552, 0.0566, 0.0577, 0.0591, 0.0602, 0.0605, 0.0609, 0.0615, 0.0627, 0.0638, 0.0644, 0.0652, 0.0661, 0.0670, 0.0678, 0.0692, 0.0706, 0.0729, 0.0757, 0.0786, 0.0805, 0.0825, 0.0846, 0.0870, 0.0897, 0.0921, 0.0947, 0.0968, 0.0997, 0.1018, 0.1027, 0.1025, 0.1018, 0.1004, 0.1000, 0.0994, 0.0980, 0.0972, 0.0960, 0.0941, 0.0927, 0.0916, 0.0902, 0.0896, 0.0890, 0.0892, 0.0896, 0.0908, 0.0919, 0.0922, 0.0937, 0.0948, 0.0957, 0.0960, 0.0961, 0.0963, 0.0965, 0.0970, 0.0983, 0.0994, 0.0997, 0.0993, 0.0984, 0.0965, 0.0951, 0.0934, 0.0916, 0.0897, 0.0870, 0.0840, 0.0813, 0.0791, 0.0766, 0.0751, 0.0730, 0.0707, 0.0683, 0.0644, 0.0616, 0.0592, 0.0562, 0.0545, 0.0531, 0.0519, 0.0504, 0.0490, 0.0468, 0.0451, 0.0432, 0.0414, 0.0403, 0.0394, 0.0386, 0.0380, 0.0370, 0.0364, 0.0367, 0.0374, 0.0385, 0.0390, 0.0390, 0.0381, 0.0380, 0.0377, 0.0381, 0.0380, 0.0377, 0.0374, 0.0376, 0.0378, 0.0380, 0.0382, 0.0385, 0.0381, 0.0377, 0.0373, 0.0367, 0.0365, 0.0358, 0.0351, 0.0342, 0.0336, 0.0334, 0.0326, 0.0322, 0.0329, 0.0327, 0.0321, 0.0310, 0.0297, 0.0293, 0.0290, 0.0283, 0.0279, 0.0272, 0.0271, 0.0271, 0.0279, 0.0282, 0.0302, 0.0325, 0.0351, 0.0375, 0.0393, 0.0406, 0.0416, 0.0422, 0.0428, 0.0430, 0.0434, 0.0443, 0.0447, 0.0457, 0.0465, 0.0479, 0.0494, 0.0514, 0.0527, 0.0539, 0.0557, 0.0571, 0.0572, 0.0563, 0.0539, 0.0504, 0.0469, 0.0441, 0.0412, 0.0385, 0.0359, 0.0334, 0.0308, 0.0282, 0.0260, 0.0232, 0.0211, 0.0196, 0.0180, 0.0169, 0.0154, 0.0137, 0.0121, 0.0105, 0.0088, 0.0067, 0.0044, 0.0027, -0.0000, -0.0024, -0.0048, -0.0066, -0.0082, -0.0111, -0.0136, -0.0158, -0.0179, -0.0201, -0.0218, -0.0235, -0.0242, -0.0245, -0.0236, -0.0231, -0.0237, -0.0237, -0.0233, -0.0229, -0.0233, -0.0239, -0.0241, -0.0244, -0.0247, -0.0251, -0.0259, -0.0270, -0.0288, -0.0295, -0.0305, -0.0311, -0.0322, -0.0327, -0.0343, -0.0352, -0.0361, -0.0358, -0.0362, -0.0365, -0.0365, -0.0358, -0.0353, -0.0348, -0.0355, -0.0365, -0.0373, -0.0373, -0.0375, -0.0365, -0.0345, -0.0327, -0.0322, -0.0327, -0.0335, -0.0337, -0.0337, -0.0348, -0.0359, -0.0361, -0.0364, -0.0371, -0.0366, -0.0361, -0.0358, -0.0353, -0.0348, -0.0345, -0.0335, -0.0320, -0.0300, -0.0281, -0.0257, -0.0233, -0.0208, -0.0179, -0.0146, -0.0104, -0.0062, -0.0031, -0.0007, 0.0023, 0.0049, 0.0077, 0.0099, 0.0125, 0.0147, 0.0177, 0.0202, 0.0232, 0.0264, 0.0291, 0.0322, 0.0346, 0.0365, 0.0386, 0.0403, 0.0415, 0.0425, 0.0428, 0.0436, 0.0447, 0.0457, 0.0465, 0.0475, 0.0494, 0.0516, 0.0534, 0.0555, 0.0574, 0.0591, 0.0616, 0.0638, 0.0655, 0.0660, 0.0657, 0.0656, 0.0650, 0.0647, 0.0637, 0.0623, 0.0613, 0.0611, 0.0611, 0.0618, 0.0633, 0.0652, 0.0665, 0.0677, 0.0690, 0.0700, 0.0712, 0.0713, 0.0710, 0.0709, 0.0695, 0.0675, 0.0655, 0.0626, 0.0598, 0.0566, 0.0533, 0.0501, 0.0470, 0.0437, 0.0405, 0.0372, 0.0342, 0.0323, 0.0305, 0.0289, 0.0267, 0.0250, 0.0229, 0.0204, 0.0183, 0.0164, 0.0152, 0.0148, 0.0141, 0.0135, 0.0139, 0.0148, 0.0160, 0.0179, 0.0193, 0.0210, 0.0234, 0.0266, 0.0291, 0.0314, 0.0337, 0.0358, 0.0375, 0.0395, 0.0412, 0.0417, 0.0424, 0.0428, 0.0428, 0.0418, 0.0415, 0.0398, 0.0373, 0.0347, 0.0328, 0.0318, 0.0301, 0.0284, 0.0258, 0.0240, 0.0214, 0.0188, 0.0168, 0.0148, 0.0137, 0.0122, 0.0109, 0.0101, 0.0093, 0.0093, 0.0096, 0.0089, 0.0104, 0.0123, 0.0141, 0.0150, 0.0159, 0.0165, 0.0174, 0.0185, 0.0205, 0.0225, 0.0247, 0.0269, 0.0290, 0.0308, 0.0331, 0.0357, 0.0379, 0.0392, 0.0405, 0.0423, 0.0440, 0.0465, 0.0484, 0.0503, 0.0506, 0.0506, 0.0504, 0.0504, 0.0513, 0.0521, 0.0537, 0.0555, 0.0577, 0.0595, 0.0619, 0.0636, 0.0654, 0.0666, 0.0674, 0.0672, 0.0667, 0.0664, 0.0660, 0.0648, 0.0645, 0.0642, 0.0645, 0.0651, 0.0653, 0.0643, 0.0629, 0.0621, 0.0607, 0.0597, 0.0585, 0.0574, 0.0553, 0.0539, 0.0525, 0.0516, 0.0508, 0.0505, 0.0500, 0.0500, 0.0490, 0.0475, 0.0464, 0.0452, 0.0440, 0.0427, 0.0416, 0.0403, 0.0393, 0.0377, 0.0360, 0.0349, 0.0343, 0.0332, 0.0325, 0.0313, 0.0305, 0.0287, 0.0263, 0.0237, 0.0211, 0.0194, 0.0182, 0.0180, 0.0184, 0.0186, 0.0186, 0.0192, 0.0198, 0.0203, 0.0200, 0.0196, 0.0177, 0.0156, 0.0130, 0.0110, 0.0087, 0.0064, 0.0049, 0.0029, 0.0011, -0.0016, -0.0042, -0.0064, -0.0080, -0.0100, -0.0118, -0.0133, -0.0147, -0.0160, -0.0170, -0.0176, -0.0188, -0.0206, -0.0209, -0.0201, -0.0200, -0.0198, -0.0191, -0.0181, -0.0170, -0.0165, -0.0161, -0.0159, -0.0157, -0.0156, -0.0156, -0.0150, -0.0137, -0.0131, -0.0128, -0.0129, -0.0124, -0.0119, -0.0106, -0.0096, -0.0086, -0.0084, -0.0075, -0.0069, -0.0070, -0.0075, -0.0092, -0.0109, -0.0124, -0.0137, -0.0140, -0.0146, -0.0146, -0.0144, -0.0144, -0.0155, -0.0160, -0.0167, -0.0174, -0.0187, -0.0198, -0.0208, -0.0217, -0.0226, -0.0229, -0.0226, -0.0219, -0.0194, -0.0173, -0.0155, -0.0143, -0.0129, -0.0120, -0.0126, -0.0137, -0.0148, -0.0153, -0.0166, -0.0182, -0.0197, -0.0216, -0.0230, -0.0244, -0.0255, -0.0270, -0.0275, -0.0280, -0.0277, -0.0276, -0.0272, -0.0272, -0.0276, -0.0284, -0.0292, -0.0300, -0.0299, -0.0305, -0.0315, -0.0316, -0.0313, -0.0320, -0.0324, -0.0326, -0.0339, -0.0354, -0.0368, -0.0376, -0.0379, -0.0381, -0.0377, -0.0373, -0.0374, -0.0380, -0.0382, -0.0394, -0.0402, -0.0415, -0.0429, -0.0439, -0.0442, -0.0455, -0.0472, -0.0487, -0.0504, -0.0515, -0.0537, -0.0555, -0.0570, -0.0583, -0.0589, -0.0601, -0.0606, -0.0614, -0.0617, -0.0621, -0.0620, -0.0614, -0.0616, -0.0618, -0.0620, -0.0622, -0.0628, -0.0623, -0.0617, -0.0602, -0.0590, -0.0584, -0.0577, -0.0572, -0.0555, -0.0536, -0.0511, -0.0488, -0.0468, -0.0448, -0.0427, -0.0401, -0.0378, -0.0358, ])
  2234. self.y[1] = array([0.0000, 0.0002, 0.0002, 0.0001, 0.0001, 0.0003, 0.0008, 0.0009, 0.0009, 0.0015, 0.0019, 0.0027, 0.0033, 0.0037, 0.0041, 0.0052, 0.0055, 0.0050, 0.0048, 0.0054, 0.0054, 0.0059, 0.0061, 0.0060, 0.0054, 0.0050, 0.0047, 0.0042, 0.0033, 0.0031, 0.0027, 0.0021, 0.0015, 0.0008, -0.0002, -0.0011, -0.0015, -0.0015, -0.0021, -0.0025, -0.0027, -0.0013, 0.0005, 0.0030, 0.0049, 0.0074, 0.0099, 0.0119, 0.0142, 0.0166, 0.0186, 0.0205, 0.0229, 0.0247, 0.0266, 0.0286, 0.0307, 0.0327, 0.0346, 0.0368, 0.0379, 0.0393, 0.0409, 0.0434, 0.0457, 0.0478, 0.0499, 0.0520, 0.0534, 0.0544, 0.0561, 0.0576, 0.0587, 0.0593, 0.0597, 0.0599, 0.0594, 0.0588, 0.0588, 0.0595, 0.0609, 0.0628, 0.0651, 0.0673, 0.0693, 0.0721, 0.0760, 0.0798, 0.0837, 0.0881, 0.0931, 0.0977, 0.1023, 0.1076, 0.1127, 0.1175, 0.1223, 0.1264, 0.1302, 0.1333, 0.1373, 0.1408, 0.1441, 0.1472, 0.1497, 0.1516, 0.1523, 0.1527, 0.1532, 0.1537, 0.1542, 0.1549, 0.1551, 0.1543, 0.1536, 0.1532, 0.1524, 0.1512, 0.1497, 0.1484, 0.1473, 0.1460, 0.1445, 0.1429, 0.1411, 0.1398, 0.1395, 0.1395, 0.1388, 0.1386, 0.1380, 0.1369, 0.1355, 0.1341, 0.1319, 0.1293, 0.1266, 0.1240, 0.1215, 0.1189, 0.1160, 0.1135, 0.1111, 0.1079, 0.1049, 0.1016, 0.0995, 0.0972, 0.0948, 0.0927, 0.0914, 0.0904, 0.0892, 0.0878, 0.0862, 0.0848, 0.0829, 0.0809, 0.0792, 0.0776, 0.0760, 0.0754, 0.0740, 0.0724, 0.0700, 0.0680, 0.0668, 0.0658, 0.0655, 0.0653, 0.0648, 0.0636, 0.0613, 0.0589, 0.0573, 0.0559, 0.0552, 0.0558, 0.0571, 0.0587, 0.0602, 0.0608, 0.0613, 0.0607, 0.0601, 0.0595, 0.0595, 0.0590, 0.0591, 0.0598, 0.0603, 0.0603, 0.0610, 0.0610, 0.0607, 0.0603, 0.0601, 0.0607, 0.0608, 0.0616, 0.0627, 0.0640, 0.0653, 0.0662, 0.0669, 0.0673, 0.0667, 0.0654, 0.0646, 0.0635, 0.0624, 0.0604, 0.0583, 0.0562, 0.0547, 0.0537, 0.0524, 0.0522, 0.0516, 0.0515, 0.0518, 0.0515, 0.0504, 0.0501, 0.0501, 0.0501, 0.0498, 0.0494, 0.0489, 0.0485, 0.0483, 0.0484, 0.0477, 0.0472, 0.0468, 0.0471, 0.0463, 0.0461, 0.0460, 0.0466, 0.0461, 0.0454, 0.0448, 0.0445, 0.0439, 0.0437, 0.0439, 0.0438, 0.0430, 0.0423, 0.0418, 0.0415, 0.0409, 0.0400, 0.0386, 0.0376, 0.0370, 0.0357, 0.0351, 0.0354, 0.0356, 0.0356, 0.0362, 0.0366, 0.0380, 0.0385, 0.0378, 0.0381, 0.0390, 0.0390, 0.0388, 0.0386, 0.0379, 0.0366, 0.0362, 0.0364, 0.0365, 0.0363, 0.0363, 0.0359, 0.0354, 0.0348, 0.0354, 0.0350, 0.0337, 0.0314, 0.0282, 0.0258, 0.0240, 0.0230, 0.0227, 0.0227, 0.0221, 0.0214, 0.0201, 0.0188, 0.0179, 0.0173, 0.0165, 0.0155, 0.0133, 0.0109, 0.0082, 0.0052, 0.0028, 0.0003, -0.0026, -0.0062, -0.0095, -0.0120, -0.0151, -0.0180, -0.0209, -0.0233, -0.0262, -0.0293, -0.0316, -0.0337, -0.0361, -0.0385, -0.0419, -0.0445, -0.0470, -0.0494, -0.0515, -0.0534, -0.0561, -0.0590, -0.0619, -0.0643, -0.0666, -0.0688, -0.0710, -0.0733, -0.0753, -0.0774, -0.0795, -0.0822, -0.0837, -0.0851, -0.0865, -0.0876, -0.0892, -0.0892, -0.0893, -0.0896, -0.0901, -0.0905, -0.0909, -0.0909, -0.0901, -0.0890, -0.0884, -0.0880, -0.0873, -0.0871, -0.0872, -0.0860, -0.0858, -0.0863, -0.0865, -0.0864, -0.0865, -0.0871, -0.0887, -0.0904, -0.0924, -0.0941, -0.0962, -0.0980, -0.0999, -0.1014, -0.1021, -0.1029, -0.1037, -0.1041, -0.1052, -0.1056, -0.1065, -0.1070, -0.1075, -0.1073, -0.1079, -0.1080, -0.1078, -0.1075, -0.1074, -0.1078, -0.1082, -0.1087, -0.1083, -0.1079, -0.1087, -0.1096, -0.1101, -0.1105, -0.1114, -0.1121, -0.1127, -0.1135, -0.1137, -0.1133, -0.1131, -0.1126, -0.1115, -0.1109, -0.1099, -0.1097, -0.1092, -0.1081, -0.1066, -0.1057, -0.1050, -0.1047, -0.1040, -0.1030, -0.1018, -0.1004, -0.0988, -0.0972, -0.0955, -0.0946, -0.0938, -0.0925, -0.0905, -0.0886, -0.0864, -0.0838, -0.0805, -0.0776, -0.0748, -0.0720, -0.0686, -0.0655, -0.0621, -0.0588, -0.0555, -0.0520, -0.0477, -0.0434, -0.0389, -0.0353, -0.0323, -0.0287, -0.0247, -0.0211, -0.0179, -0.0149, -0.0123, -0.0094, -0.0065, -0.0042, -0.0020, -0.0001, 0.0016, 0.0028, 0.0044, 0.0072, 0.0097, 0.0109, 0.0120, 0.0127, 0.0135, 0.0140, 0.0151, 0.0147, 0.0158, 0.0161, 0.0162, 0.0154, 0.0154, 0.0162, 0.0172, 0.0186, 0.0204, 0.0214, 0.0213, 0.0210, 0.0211, 0.0208, 0.0204, 0.0202, 0.0201, 0.0187, 0.0179, 0.0167, 0.0164, 0.0167, 0.0170, 0.0180, 0.0188, 0.0199, 0.0204, 0.0213, 0.0209, 0.0203, 0.0198, 0.0186, 0.0183, 0.0181, 0.0168, 0.0163, 0.0161, 0.0160, 0.0154, 0.0157, 0.0162, 0.0166, 0.0171, 0.0171, 0.0179, 0.0186, 0.0196, 0.0205, 0.0229, 0.0254, 0.0280, 0.0305, 0.0327, 0.0343, 0.0354, 0.0361, 0.0367, 0.0363, 0.0361, 0.0355, 0.0353, 0.0351, 0.0342, 0.0335, 0.0328, 0.0326, 0.0324, 0.0325, 0.0330, 0.0333, 0.0333, 0.0327, 0.0326, 0.0328, 0.0326, 0.0327, 0.0329, 0.0328, 0.0332, 0.0331, 0.0325, 0.0315, 0.0311, 0.0305, 0.0302, 0.0300, 0.0289, 0.0275, 0.0260, 0.0246, 0.0234, 0.0230, 0.0216, 0.0211, 0.0205, 0.0191, 0.0172, 0.0147, 0.0131, 0.0112, 0.0094, 0.0082, 0.0062, 0.0045, 0.0041, 0.0035, 0.0037, 0.0035, 0.0033, 0.0033, 0.0037, 0.0034, 0.0036, 0.0033, 0.0028, 0.0021, 0.0015, 0.0008, 0.0004, 0.0002, -0.0002, -0.0008, -0.0001, 0.0007, 0.0015, 0.0020, 0.0020, 0.0028, 0.0037, 0.0045, 0.0052, 0.0060, 0.0063, 0.0074, 0.0083, 0.0092, 0.0095, 0.0095, 0.0096, 0.0093, 0.0100, 0.0097, 0.0093, 0.0084, 0.0071, 0.0050, 0.0044, 0.0031, 0.0017, 0.0005, -0.0004, -0.0018, -0.0026, -0.0026, -0.0020, -0.0018, -0.0010, 0.0003, 0.0024, 0.0046, 0.0066, 0.0090, 0.0110, 0.0126, 0.0143, 0.0157, 0.0167, 0.0171, 0.0179, 0.0195, 0.0223, 0.0254, 0.0287, 0.0320, 0.0359, 0.0400, 0.0439, 0.0470, 0.0491, 0.0515, 0.0535, 0.0563, 0.0585, 0.0605, 0.0620, 0.0635, 0.0654, 0.0670, 0.0687, 0.0701, 0.0728, 0.0746, 0.0766, 0.0785, 0.0800, 0.0812, 0.0817, 0.0822, 0.0820, 0.0814, 0.0814, 0.0817, 0.0820, 0.0823, 0.0826, 0.0826, 0.0831, 0.0833, 0.0846, 0.0851, 0.0853, 0.0855, 0.0858, 0.0859, 0.0864, 0.0868, 0.0876, 0.0884, 0.0890, 0.0894, 0.0895, 0.0888, 0.0882, 0.0883, 0.0884, 0.0892, 0.0903, 0.0915, 0.0917, 0.0916, 0.0908, 0.0900, 0.0891, 0.0881, 0.0869, 0.0849, 0.0825, 0.0805, 0.0789, 0.0773, 0.0761, 0.0750, 0.0738, 0.0729, 0.0717, 0.0712, 0.0714, 0.0718, 0.0727, 0.0737, 0.0735, 0.0733, 0.0723, 0.0706, 0.0683, 0.0666, 0.0651, 0.0636, 0.0620, 0.0599, 0.0580, 0.0563, 0.0551, 0.0540, 0.0526, 0.0516, 0.0503, 0.0492, 0.0482, 0.0469, 0.0468, 0.0471, 0.0471, 0.0471, 0.0483, 0.0490, 0.0486, 0.0489, 0.0495, 0.0510, 0.0523, 0.0538, 0.0547, 0.0552, 0.0556, 0.0569, 0.0573, 0.0578, 0.0583, 0.0579, 0.0579, 0.0580, 0.0583, 0.0585, 0.0594, 0.0601, 0.0598, 0.0594, 0.0597, 0.0602, 0.0611, 0.0612, 0.0618, 0.0620, 0.0626, 0.0635, 0.0637, 0.0645, 0.0649, 0.0657, 0.0667, 0.0680, 0.0688, 0.0693, 0.0698, 0.0702, 0.0699, 0.0703, 0.0702, 0.0705, 0.0716, 0.0725, 0.0730, 0.0736, 0.0736, 0.0734, 0.0735, 0.0739, 0.0736, 0.0728, 0.0720, 0.0717, 0.0715, 0.0709, 0.0704, 0.0707, 0.0704, 0.0695, 0.0690, 0.0687, 0.0667, 0.0648, 0.0629, 0.0620, 0.0620, 0.0618, 0.0621, 0.0623, 0.0633, 0.0637, 0.0637, 0.0638, 0.0639, 0.0639, 0.0642, 0.0650, 0.0653, 0.0657, 0.0661, ])
  2235. self.y[2] = array([0.0000, 0.0001, -0.0002, -0.0003, 0.0004, 0.0014, 0.0021, 0.0025, 0.0025, 0.0021, 0.0018, 0.0022, 0.0016, 0.0018, 0.0018, 0.0021, 0.0027, 0.0034, 0.0046, 0.0060, 0.0076, 0.0080, 0.0084, 0.0090, 0.0100, 0.0104, 0.0098, 0.0097, 0.0100, 0.0100, 0.0105, 0.0117, 0.0124, 0.0128, 0.0133, 0.0133, 0.0133, 0.0132, 0.0132, 0.0136, 0.0144, 0.0161, 0.0179, 0.0196, 0.0222, 0.0251, 0.0265, 0.0279, 0.0287, 0.0291, 0.0297, 0.0305, 0.0316, 0.0328, 0.0340, 0.0361, 0.0382, 0.0408, 0.0425, 0.0442, 0.0460, 0.0474, 0.0489, 0.0502, 0.0512, 0.0517, 0.0526, 0.0525, 0.0525, 0.0521, 0.0508, 0.0498, 0.0487, 0.0478, 0.0472, 0.0461, 0.0446, 0.0434, 0.0415, 0.0399, 0.0388, 0.0374, 0.0360, 0.0351, 0.0339, 0.0320, 0.0308, 0.0294, 0.0286, 0.0278, 0.0255, 0.0224, 0.0194, 0.0170, 0.0147, 0.0131, 0.0123, 0.0121, 0.0110, 0.0106, 0.0100, 0.0090, 0.0089, 0.0093, 0.0100, 0.0111, 0.0132, 0.0159, 0.0179, 0.0195, 0.0207, 0.0220, 0.0229, 0.0242, 0.0261, 0.0279, 0.0290, 0.0300, 0.0303, 0.0309, 0.0316, 0.0328, 0.0335, 0.0332, 0.0323, 0.0314, 0.0300, 0.0283, 0.0268, 0.0248, 0.0225, 0.0201, 0.0172, 0.0143, 0.0119, 0.0104, 0.0088, 0.0071, 0.0059, 0.0051, 0.0038, 0.0027, 0.0017, 0.0009, 0.0007, -0.0002, -0.0008, -0.0020, -0.0034, -0.0055, -0.0069, -0.0082, -0.0086, -0.0085, -0.0089, -0.0088, -0.0092, -0.0099, -0.0110, -0.0121, -0.0137, -0.0150, -0.0157, -0.0160, -0.0155, -0.0149, -0.0140, -0.0139, -0.0134, -0.0135, -0.0136, -0.0138, -0.0145, -0.0158, -0.0155, -0.0155, -0.0153, -0.0154, -0.0159, -0.0157, -0.0152, -0.0143, -0.0140, -0.0138, -0.0142, -0.0148, -0.0154, -0.0154, -0.0147, -0.0136, -0.0122, -0.0107, -0.0102, -0.0091, -0.0082, -0.0063, -0.0042, -0.0028, -0.0007, 0.0015, 0.0038, 0.0054, 0.0075, 0.0093, 0.0110, 0.0137, 0.0156, 0.0170, 0.0184, 0.0195, 0.0199, 0.0206, 0.0209, 0.0207, 0.0201, 0.0198, 0.0196, 0.0191, 0.0192, 0.0193, 0.0196, 0.0201, 0.0206, 0.0216, 0.0229, 0.0253, 0.0278, 0.0304, 0.0321, 0.0341, 0.0356, 0.0367, 0.0372, 0.0379, 0.0383, 0.0394, 0.0396, 0.0397, 0.0391, 0.0382, 0.0362, 0.0346, 0.0334, 0.0325, 0.0325, 0.0325, 0.0323, 0.0316, 0.0321, 0.0328, 0.0338, 0.0352, 0.0366, 0.0378, 0.0389, 0.0401, 0.0403, 0.0406, 0.0419, 0.0425, 0.0424, 0.0421, 0.0415, 0.0407, 0.0406, 0.0411, 0.0413, 0.0419, 0.0420, 0.0416, 0.0415, 0.0402, 0.0388, 0.0377, 0.0368, 0.0357, 0.0353, 0.0351, 0.0348, 0.0347, 0.0340, 0.0326, 0.0314, 0.0306, 0.0295, 0.0288, 0.0275, 0.0257, 0.0234, 0.0209, 0.0188, 0.0169, 0.0151, 0.0132, 0.0114, 0.0093, 0.0076, 0.0055, 0.0026, -0.0008, -0.0038, -0.0069, -0.0097, -0.0118, -0.0132, -0.0148, -0.0164, -0.0183, -0.0202, -0.0213, -0.0228, -0.0253, -0.0280, -0.0299, -0.0315, -0.0329, -0.0343, -0.0349, -0.0364, -0.0373, -0.0383, -0.0391, -0.0400, -0.0401, -0.0404, -0.0414, -0.0422, -0.0425, -0.0425, -0.0420, -0.0419, -0.0412, -0.0412, -0.0408, -0.0400, -0.0396, -0.0397, -0.0400, -0.0401, -0.0405, -0.0412, -0.0410, -0.0401, -0.0389, -0.0377, -0.0372, -0.0369, -0.0364, -0.0356, -0.0347, -0.0345, -0.0346, -0.0353, -0.0359, -0.0367, -0.0376, -0.0385, -0.0387, -0.0391, -0.0401, -0.0400, -0.0398, -0.0395, -0.0391, -0.0385, -0.0383, -0.0377, -0.0374, -0.0371, -0.0361, -0.0357, -0.0355, -0.0342, -0.0324, -0.0303, -0.0280, -0.0257, -0.0232, -0.0205, -0.0179, -0.0147, -0.0113, -0.0091, -0.0059, -0.0019, 0.0024, 0.0071, 0.0120, 0.0177, 0.0219, 0.0248, 0.0273, 0.0290, 0.0302, 0.0318, 0.0325, 0.0331, 0.0334, 0.0344, 0.0360, 0.0387, 0.0414, 0.0437, 0.0456, 0.0471, 0.0491, 0.0511, 0.0526, 0.0536, 0.0550, 0.0570, 0.0597, 0.0622, 0.0647, 0.0672, 0.0696, 0.0710, 0.0725, 0.0748, 0.0777, 0.0806, 0.0832, 0.0854, 0.0874, 0.0893, 0.0906, 0.0914, 0.0927, 0.0937, 0.0944, 0.0959, 0.0970, 0.0970, 0.0967, 0.0955, 0.0945, 0.0934, 0.0927, 0.0921, 0.0916, 0.0909, 0.0901, 0.0893, 0.0884, 0.0881, 0.0877, 0.0870, 0.0870, 0.0871, 0.0872, 0.0867, 0.0860, 0.0845, 0.0821, 0.0800, 0.0775, 0.0752, 0.0728, 0.0703, 0.0676, 0.0655, 0.0637, 0.0618, 0.0593, 0.0558, 0.0518, 0.0483, 0.0445, 0.0414, 0.0384, 0.0365, 0.0349, 0.0338, 0.0330, 0.0323, 0.0320, 0.0319, 0.0320, 0.0323, 0.0324, 0.0324, 0.0322, 0.0327, 0.0328, 0.0320, 0.0314, 0.0295, 0.0279, 0.0258, 0.0243, 0.0229, 0.0211, 0.0195, 0.0176, 0.0163, 0.0151, 0.0136, 0.0123, 0.0108, 0.0088, 0.0064, 0.0041, 0.0021, 0.0002, -0.0020, -0.0051, -0.0078, -0.0106, -0.0137, -0.0169, -0.0195, -0.0218, -0.0235, -0.0253, -0.0269, -0.0281, -0.0290, -0.0299, -0.0315, -0.0328, -0.0336, -0.0344, -0.0352, -0.0354, -0.0353, -0.0358, -0.0355, -0.0347, -0.0344, -0.0338, -0.0335, -0.0324, -0.0319, -0.0319, -0.0322, -0.0334, -0.0351, -0.0365, -0.0381, -0.0401, -0.0424, -0.0449, -0.0472, -0.0491, -0.0506, -0.0522, -0.0533, -0.0544, -0.0553, -0.0553, -0.0557, -0.0558, -0.0563, -0.0573, -0.0587, -0.0602, -0.0617, -0.0636, -0.0666, -0.0699, -0.0734, -0.0760, -0.0791, -0.0817, -0.0840, -0.0857, -0.0869, -0.0881, -0.0882, -0.0884, -0.0888, -0.0893, -0.0889, -0.0876, -0.0869, -0.0858, -0.0851, -0.0836, -0.0825, -0.0813, -0.0807, -0.0800, -0.0802, -0.0805, -0.0807, -0.0813, -0.0818, -0.0819, -0.0822, -0.0829, -0.0825, -0.0829, -0.0838, -0.0840, -0.0836, -0.0824, -0.0810, -0.0799, -0.0789, -0.0785, -0.0785, -0.0774, -0.0770, -0.0761, -0.0755, -0.0752, -0.0746, -0.0744, -0.0741, -0.0729, -0.0718, -0.0710, -0.0694, -0.0673, -0.0649, -0.0627, -0.0606, -0.0592, -0.0582, -0.0571, -0.0555, -0.0539, -0.0527, -0.0516, -0.0508, -0.0503, -0.0508, -0.0518, -0.0527, -0.0536, -0.0544, -0.0545, -0.0541, -0.0547, -0.0553, -0.0561, -0.0573, -0.0587, -0.0600, -0.0608, -0.0616, -0.0624, -0.0626, -0.0631, -0.0638, -0.0630, -0.0621, -0.0612, -0.0603, -0.0594, -0.0577, -0.0564, -0.0551, -0.0535, -0.0511, -0.0485, -0.0455, -0.0424, -0.0396, -0.0374, -0.0360, -0.0351, -0.0343, -0.0336, -0.0324, -0.0311, -0.0296, -0.0276, -0.0259, -0.0240, -0.0223, -0.0212, -0.0205, -0.0201, -0.0203, -0.0214, -0.0232, -0.0261, -0.0285, -0.0304, -0.0320, -0.0342, -0.0365, -0.0381, -0.0394, -0.0416, -0.0434, -0.0443, -0.0450, -0.0465, -0.0482, -0.0490, -0.0497, -0.0508, -0.0506, -0.0509, -0.0514, -0.0522, -0.0525, -0.0532, -0.0547, -0.0553, -0.0575, -0.0597, -0.0629, -0.0666, -0.0702, -0.0737, -0.0773, -0.0805, -0.0835, -0.0847, -0.0875, -0.0893, -0.0913, -0.0929, -0.0935, -0.0943, -0.0949, -0.0957, -0.0962, -0.0977, -0.0994, -0.1011, -0.1031, -0.1057, -0.1086, -0.1114, -0.1141, -0.1158, -0.1174, -0.1197, -0.1219, -0.1234, -0.1250, -0.1267, -0.1273, -0.1270, -0.1268, -0.1263, -0.1257, -0.1240, -0.1220, -0.1210, -0.1197, -0.1197, -0.1188, -0.1171, -0.1147, -0.1122, -0.1101, -0.1077, -0.1052, -0.1031, -0.1020, -0.1008, -0.0983, -0.0960, -0.0934, -0.0908, -0.0889, -0.0865, -0.0836, -0.0799, -0.0767, -0.0740, -0.0716, -0.0692, -0.0663, -0.0628, -0.0596, -0.0567, -0.0535, -0.0508, -0.0472, -0.0435, -0.0394, -0.0347, -0.0301, -0.0253, -0.0210, -0.0165, -0.0126, -0.0095, -0.0072, -0.0046, -0.0020, -0.0001, 0.0023, 0.0050, 0.0077, 0.0102, 0.0135, 0.0175, 0.0212, 0.0245, 0.0280, 0.0309, 0.0336, 0.0365, 0.0397, 0.0433, 0.0474, 0.0511, 0.0546, 0.0574, 0.0602, 0.0634, 0.0663, 0.0695, 0.0729, 0.0752, 0.0762, 0.0773, 0.0781, 0.0790, 0.0806, 0.0836, 0.0857, 0.0879, 0.0896, 0.0920, 0.0949, 0.0975, 0.1002, ])
  2236. # COMPOSITE types:
  2237. # MassSpringForce: Line(horizontal), Spring, Rectangle, Arrow/Line(w/arrow)
  2238. # must be easy to find the tip of the arrow
  2239. # Maybe extra dict: self.name['mass'] = Rectangle object - YES!
  2240. def _test1():
  2241. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2242. l1 = Line((0,0), (1,1))
  2243. l1.draw()
  2244. eval(input(': '))
  2245. c1 = Circle((5,2), 1)
  2246. c2 = Circle((6,2), 1)
  2247. w1 = Wheel((7,2), 1)
  2248. c1.draw()
  2249. c2.draw()
  2250. w1.draw()
  2251. hardcopy()
  2252. display() # show the plot
  2253. def _test2():
  2254. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2255. l1 = Line((0,0), (1,1))
  2256. l1.draw()
  2257. eval(input(': '))
  2258. c1 = Circle((5,2), 1)
  2259. c2 = Circle((6,2), 1)
  2260. w1 = Wheel((7,2), 1)
  2261. filled_curves(True)
  2262. set_linecolor('blue')
  2263. c1.draw()
  2264. set_linecolor('aqua')
  2265. c2.draw()
  2266. filled_curves(False)
  2267. set_linecolor('red')
  2268. w1.draw()
  2269. hardcopy()
  2270. display() # show the plot
  2271. def _test3():
  2272. """Test example from the book."""
  2273. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2274. l1 = Line(start=(0,0), stop=(1,1)) # define line
  2275. l1.draw() # make plot data
  2276. r1 = Rectangle(lower_left_corner=(0,1), width=3, height=5)
  2277. r1.draw()
  2278. Circle(center=(5,7), radius=1).draw()
  2279. Wheel(center=(6,2), radius=2, inner_radius=0.5, nlines=7).draw()
  2280. hardcopy()
  2281. display()
  2282. def _test4():
  2283. """Second example from the book."""
  2284. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2285. r1 = Rectangle(lower_left_corner=(0,1), width=3, height=5)
  2286. c1 = Circle(center=(5,7), radius=1)
  2287. w1 = Wheel(center=(6,2), radius=2, inner_radius=0.5, nlines=7)
  2288. c2 = Circle(center=(7,7), radius=1)
  2289. filled_curves(True)
  2290. c1.draw()
  2291. set_linecolor('blue')
  2292. r1.draw()
  2293. set_linecolor('aqua')
  2294. c2.draw()
  2295. # Add thick aqua line around rectangle:
  2296. filled_curves(False)
  2297. set_linewidth(4)
  2298. r1.draw()
  2299. set_linecolor('red')
  2300. # Draw wheel with thick lines:
  2301. w1.draw()
  2302. hardcopy('tmp_colors')
  2303. display()
  2304. def _test5():
  2305. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2306. c = 6. # center point of box
  2307. w = 2. # size of box
  2308. L = 3
  2309. r1 = Rectangle((c-w/2, c-w/2), w, w)
  2310. l1 = Line((c,c-w/2), (c,c-w/2-L))
  2311. linecolor('blue')
  2312. filled_curves(True)
  2313. r1.draw()
  2314. linecolor('aqua')
  2315. filled_curves(False)
  2316. l1.draw()
  2317. hardcopy()
  2318. display() # show the plot
  2319. def rolling_wheel(total_rotation_angle):
  2320. """Animation of a rotating wheel."""
  2321. set_coordinate_system(xmin=0, xmax=10, ymin=0, ymax=10)
  2322. import time
  2323. center = (6,2)
  2324. radius = 2.0
  2325. angle = 2.0
  2326. pngfiles = []
  2327. w1 = Wheel(center=center, radius=radius, inner_radius=0.5, nlines=7)
  2328. for i in range(int(total_rotation_angle/angle)):
  2329. w1.draw()
  2330. print('BIG PROBLEM WITH ANIMATE!!!')
  2331. display()
  2332. filename = 'tmp_%03d' % i
  2333. pngfiles.append(filename + '.png')
  2334. hardcopy(filename)
  2335. time.sleep(0.3) # pause
  2336. L = radius*angle*pi/180 # translation = arc length
  2337. w1.rotate(angle, center)
  2338. w1.translate((-L, 0))
  2339. center = (center[0] - L, center[1])
  2340. erase()
  2341. cmd = 'convert -delay 50 -loop 1000 %s tmp_movie.gif' \
  2342. % (' '.join(pngfiles))
  2343. print('converting PNG files to animated GIF:\n', cmd)
  2344. import subprocess
  2345. failure, output = subprocess.getstatusoutput(cmd)
  2346. if failure: print('Could not run', cmd)
  2347. if __name__ == '__main__':
  2348. #rolling_wheel(40)
  2349. #_test1()
  2350. #_test3()
  2351. funcs = [
  2352. #test_Axis,
  2353. test_inclined_plane,
  2354. ]
  2355. for func in funcs:
  2356. func()
  2357. input('Type Return: ')