1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31 """ Implementation of the RECAP algorithm from Lewell et al. JCICS *38* 511-522 (1998)
32
33 The published algorithm is implemented more or less without
34 modification. The results are returned as a hierarchy of nodes instead
35 of just as a set of fragments. The hope is that this will allow a bit
36 more flexibility in working with the results.
37
38 For example:
39 >>> from rdkit import Chem
40 >>> from rdkit.Chem import Recap
41 >>> m = Chem.MolFromSmiles('C1CC1Oc1ccccc1-c1ncc(OC)cc1')
42 >>> res = Recap.RecapDecompose(m)
43 >>> res
44 <...Chem.Recap.RecapHierarchyNode object at ...>
45 >>> sorted(res.children.keys())
46 ['[*]C1CC1', '[*]c1ccc(OC)cn1', '[*]c1ccccc1-c1...cc(OC)c...1', '[*]c1ccccc1OC1CC1']
47 >>> sorted(res.GetAllChildren().keys())
48 ['[*]C1CC1', '[*]c1ccc(OC)cn1', '[*]c1ccccc1-c1...cc(OC)c...1', '[*]c1ccccc1OC1CC1', '[*]c1ccccc1[*]']
49
50
51 To get the standard set of RECAP results, use GetLeaves():
52 >>> leaves=res.GetLeaves()
53 >>> sorted(leaves.keys())
54 ['[*]C1CC1', '[*]c1ccc(OC)cn1', '[*]c1ccccc1[*]']
55 >>> leaf = leaves['[*]C1CC1']
56 >>> leaf.mol
57 <...Chem.rdchem.Mol object at ...>
58
59
60 """
61 import sys
62 import weakref
63 from rdkit import Chem
64 from rdkit.Chem import rdChemReactions as Reactions
65 from rdkit.six import iterkeys, iteritems, next
66
67
68 reactionDefs = (
69 "[#7;+0;D2,D3:1]!@C(!@=O)!@[#7;+0;D2,D3:2]>>[*][#7:1].[#7:2][*]",
70 "[C;!$(C([#7])[#7]):1](=!@[O:2])!@[#7;+0;!D1:3]>>[*][C:1]=[O:2].[*][#7:3]",
71 "[C:1](=!@[O:2])!@[O;+0:3]>>[*][C:1]=[O:2].[O:3][*]",
72 "[N;!D1;+0;!$(N-C=[#7,#8,#15,#16])](-!@[*:1])-!@[*:2]>>[*][*:1].[*:2][*]",
73
74
75
76 "[#7;R;D3;+0:1]-!@[*:2]>>[*][#7:1].[*:2][*]",
77 "[#6:1]-!@[O;+0]-!@[#6:2]>>[#6:1][*].[*][#6:2]",
78 "[C:1]=!@[C:2]>>[C:1][*].[*][C:2]",
79 "[n;+0:1]-!@[C:2]>>[n:1][*].[C:2][*]",
80 "[O:3]=[C:4]-@[N;+0:1]-!@[C:2]>>[O:3]=[C:4]-[N:1][*].[C:2][*]",
81 "[c:1]-!@[c:2]>>[c:1][*].[*][c:2]",
82 "[n;+0:1]-!@[c:2]>>[n:1][*].[*][c:2]",
83 "[#7;+0;D2,D3:1]-!@[S:2](=[O:3])=[O:4]>>[#7:1][*].[*][S:2](=[O:3])=[O:4]",
84 )
85
86 reactions = tuple([Reactions.ReactionFromSmarts(x) for x in reactionDefs])
87
88
90 """ This class is used to hold the Recap hiearchy
91 """
92 mol = None
93 children = None
94 parents = None
95 smiles = None
96
101
103 " returns a dictionary, keyed by SMILES, of children "
104 res = {}
105 for smi, child in iteritems(self.children):
106 res[smi] = child
107 child._gacRecurse(res, terminalOnly=False)
108 return res
109
111 " returns a dictionary, keyed by SMILES, of leaf (terminal) nodes "
112 res = {}
113 for smi, child in iteritems(self.children):
114 if not len(child.children):
115 res[smi] = child
116 else:
117 child._gacRecurse(res, terminalOnly=True)
118 return res
119
121 """ returns all the nodes in the hierarchy tree that contain this
122 node as a child
123 """
124 if not self.parents:
125 res = [self]
126 else:
127 res = []
128 for p in self.parents.values():
129 for uP in p.getUltimateParents():
130 if uP not in res:
131 res.append(uP)
132 return res
133
135 for smi, child in iteritems(self.children):
136 if not terminalOnly or not len(child.children):
137 res[smi] = child
138 child._gacRecurse(res, terminalOnly=terminalOnly)
139
144
145
146 -def RecapDecompose(mol, allNodes=None, minFragmentSize=0, onlyUseReactions=None):
147 """ returns the recap decomposition for a molecule """
148 mSmi = Chem.MolToSmiles(mol, 1)
149
150 if allNodes is None:
151 allNodes = {}
152 if mSmi in allNodes:
153 return allNodes[mSmi]
154
155 res = RecapHierarchyNode(mol)
156 res.smiles = mSmi
157 activePool = {mSmi: res}
158 allNodes[mSmi] = res
159 while activePool:
160 nSmi = next(iterkeys(activePool))
161 node = activePool.pop(nSmi)
162 if not node.mol:
163 continue
164 for rxnIdx, reaction in enumerate(reactions):
165 if onlyUseReactions and rxnIdx not in onlyUseReactions:
166 continue
167
168
169 ps = reaction.RunReactants((node.mol, ))
170
171 if ps:
172 for prodSeq in ps:
173 seqOk = True
174
175
176 tSeq = [(prod.GetNumAtoms(onlyExplicit=True), idx) for idx, prod in enumerate(prodSeq)]
177 tSeq.sort()
178 ts = [(x, prodSeq[y]) for x, y in tSeq]
179 prodSeq = ts
180 for nats, prod in prodSeq:
181 try:
182 Chem.SanitizeMol(prod)
183 except Exception:
184 continue
185 pSmi = Chem.MolToSmiles(prod, 1)
186 if minFragmentSize > 0:
187 nDummies = pSmi.count('*')
188 if nats - nDummies < minFragmentSize:
189 seqOk = False
190 break
191
192
193 elif pSmi.replace('[*]', '').replace('()', '') in ('', 'C', 'CC', 'CCC'):
194 seqOk = False
195 break
196 prod.pSmi = pSmi
197 if seqOk:
198 for nats, prod in prodSeq:
199 pSmi = prod.pSmi
200
201 if not pSmi in allNodes:
202 pNode = RecapHierarchyNode(prod)
203 pNode.smiles = pSmi
204 pNode.parents[nSmi] = weakref.proxy(node)
205 node.children[pSmi] = pNode
206 activePool[pSmi] = pNode
207 allNodes[pSmi] = pNode
208 else:
209 pNode = allNodes[pSmi]
210 pNode.parents[nSmi] = weakref.proxy(node)
211 node.children[pSmi] = pNode
212
213 return res
214
215
216
217
218 if __name__ == '__main__':
219 import unittest
220
222
224 m = Chem.MolFromSmiles('C1CC1Oc1ccccc1-c1ncc(OC)cc1')
225 res = RecapDecompose(m)
226 self.assertTrue(res)
227 self.assertTrue(len(res.children.keys()) == 4)
228 self.assertTrue(len(res.GetAllChildren().keys()) == 5)
229 self.assertTrue(len(res.GetLeaves().keys()) == 3)
230
232 m = Chem.MolFromSmiles('CCCOCCC')
233 res = RecapDecompose(m)
234 self.assertTrue(res)
235 self.assertTrue(res.children == {})
236
238 allNodes = {}
239 m = Chem.MolFromSmiles('c1ccccc1-c1ncccc1')
240 res = RecapDecompose(m, allNodes=allNodes)
241 self.assertTrue(res)
242 self.assertTrue(len(res.children.keys()) == 2)
243 self.assertTrue(len(allNodes.keys()) == 3)
244
245 m = Chem.MolFromSmiles('COc1ccccc1-c1ncccc1')
246 res = RecapDecompose(m, allNodes=allNodes)
247 self.assertTrue(res)
248 self.assertTrue(len(res.children.keys()) == 2)
249
250 self.assertTrue(len(allNodes.keys()) == 5)
251 self.assertTrue('[*]c1ccccc1OC' in allNodes)
252 self.assertTrue('[*]c1ccccc1' in allNodes)
253
254 m = Chem.MolFromSmiles('C1CC1Oc1ccccc1-c1ncccc1')
255 res = RecapDecompose(m, allNodes=allNodes)
256 self.assertTrue(res)
257 self.assertTrue(len(res.children.keys()) == 4)
258 self.assertTrue(len(allNodes.keys()) == 10)
259
261 m = Chem.MolFromSmiles('c1ccccc1OC(Oc1ccccc1)Oc1ccccc1')
262 res = RecapDecompose(m)
263 self.assertTrue(res)
264 self.assertTrue(len(res.GetLeaves()) == 2)
265 ks = res.GetLeaves().keys()
266 self.assertFalse('[*]C([*])[*]' in ks)
267 self.assertTrue('[*]c1ccccc1' in ks)
268 self.assertTrue('[*]C([*])Oc1ccccc1' in ks)
269
271 m = Chem.MolFromSmiles('C1CCCCN1CCCC')
272 res = RecapDecompose(m)
273 self.assertTrue(res)
274 self.assertTrue(len(res.GetLeaves()) == 2)
275 ks = res.GetLeaves().keys()
276 self.assertTrue('[*]N1CCCCC1' in ks)
277 self.assertTrue('[*]CCCC' in ks)
278
280 m = Chem.MolFromSmiles('CCCOCCC')
281 res = RecapDecompose(m)
282 self.assertTrue(res)
283 self.assertTrue(res.children == {})
284 res = RecapDecompose(m, minFragmentSize=3)
285 self.assertTrue(res)
286 self.assertTrue(len(res.GetLeaves()) == 1)
287 ks = res.GetLeaves().keys()
288 self.assertTrue('[*]CCC' in ks)
289
290 m = Chem.MolFromSmiles('CCCOCC')
291 res = RecapDecompose(m, minFragmentSize=3)
292 self.assertTrue(res)
293 self.assertTrue(res.children == {})
294
295 m = Chem.MolFromSmiles('CCCOCCOC')
296 res = RecapDecompose(m, minFragmentSize=2)
297 self.assertTrue(res)
298 self.assertTrue(len(res.GetLeaves()) == 2)
299 ks = res.GetLeaves().keys()
300 self.assertTrue('[*]CCC' in ks)
301 ks = res.GetLeaves().keys()
302 self.assertTrue('[*]CCOC' in ks)
303
305 m = Chem.MolFromSmiles('C1CC1C(=O)NC1OC1')
306 res = RecapDecompose(m, onlyUseReactions=[1])
307 self.assertTrue(res)
308 self.assertTrue(len(res.GetLeaves()) == 2)
309 ks = res.GetLeaves().keys()
310 self.assertTrue('[*]C(=O)C1CC1' in ks)
311 self.assertTrue('[*]NC1CO1' in ks)
312
313 m = Chem.MolFromSmiles('C1CC1C(=O)N(C)C1OC1')
314 res = RecapDecompose(m, onlyUseReactions=[1])
315 self.assertTrue(res)
316 self.assertTrue(len(res.GetLeaves()) == 2)
317 ks = res.GetLeaves().keys()
318 self.assertTrue('[*]C(=O)C1CC1' in ks)
319 self.assertTrue('[*]N(C)C1CO1' in ks)
320
321 m = Chem.MolFromSmiles('C1CC1C(=O)n1cccc1')
322 res = RecapDecompose(m, onlyUseReactions=[1])
323 self.assertTrue(res)
324 self.assertTrue(len(res.GetLeaves()) == 2)
325 ks = res.GetLeaves().keys()
326 self.assertTrue('[*]C(=O)C1CC1' in ks)
327 self.assertTrue('[*]n1cccc1' in ks)
328
329 m = Chem.MolFromSmiles('C1CC1C(=O)CC1OC1')
330 res = RecapDecompose(m, onlyUseReactions=[1])
331 self.assertTrue(res)
332 self.assertTrue(len(res.GetLeaves()) == 0)
333
334 m = Chem.MolFromSmiles('C1CCC(=O)NC1')
335 res = RecapDecompose(m, onlyUseReactions=[1])
336 self.assertTrue(res)
337 self.assertTrue(len(res.GetLeaves()) == 0)
338
339 m = Chem.MolFromSmiles('CC(=O)NC')
340 res = RecapDecompose(m, onlyUseReactions=[1])
341 self.assertTrue(res)
342 self.assertTrue(len(res.GetLeaves()) == 2)
343 ks = res.GetLeaves().keys()
344
345 m = Chem.MolFromSmiles('CC(=O)N')
346 res = RecapDecompose(m, onlyUseReactions=[1])
347 self.assertTrue(res)
348 self.assertTrue(len(res.GetLeaves()) == 0)
349
350 m = Chem.MolFromSmiles('C(=O)NCCNC(=O)CC')
351 res = RecapDecompose(m, onlyUseReactions=[1])
352 self.assertTrue(res)
353 self.assertTrue(len(res.children) == 4)
354 self.assertTrue(len(res.GetLeaves()) == 3)
355
357 m = Chem.MolFromSmiles('C1CC1C(=O)OC1OC1')
358 res = RecapDecompose(m, onlyUseReactions=[2])
359 self.assertTrue(res)
360 self.assertTrue(len(res.GetLeaves()) == 2)
361 ks = res.GetLeaves().keys()
362 self.assertTrue('[*]C(=O)C1CC1' in ks)
363 self.assertTrue('[*]OC1CO1' in ks)
364
365 m = Chem.MolFromSmiles('C1CC1C(=O)CC1OC1')
366 res = RecapDecompose(m, onlyUseReactions=[2])
367 self.assertTrue(res)
368 self.assertTrue(len(res.GetLeaves()) == 0)
369
370 m = Chem.MolFromSmiles('C1CCC(=O)OC1')
371 res = RecapDecompose(m, onlyUseReactions=[2])
372 self.assertTrue(res)
373 self.assertTrue(len(res.GetLeaves()) == 0)
374
376 m = Chem.MolFromSmiles('C1CC1NC(=O)NC1OC1')
377 res = RecapDecompose(m, onlyUseReactions=[0])
378 self.assertTrue(res)
379 self.assertTrue(len(res.GetLeaves()) == 2)
380 ks = res.GetLeaves().keys()
381 self.assertTrue('[*]NC1CC1' in ks)
382 self.assertTrue('[*]NC1CO1' in ks)
383
384 m = Chem.MolFromSmiles('C1CC1NC(=O)N(C)C1OC1')
385 res = RecapDecompose(m, onlyUseReactions=[0])
386 self.assertTrue(res)
387 self.assertTrue(len(res.GetLeaves()) == 2)
388 ks = res.GetLeaves().keys()
389 self.assertTrue('[*]NC1CC1' in ks)
390 self.assertTrue('[*]N(C)C1CO1' in ks)
391
392 m = Chem.MolFromSmiles('C1CCNC(=O)NC1C')
393 res = RecapDecompose(m, onlyUseReactions=[0])
394 self.assertTrue(res)
395 self.assertTrue(len(res.GetLeaves()) == 0)
396
397 m = Chem.MolFromSmiles('c1cccn1C(=O)NC1OC1')
398 res = RecapDecompose(m, onlyUseReactions=[0])
399 self.assertTrue(res)
400 self.assertTrue(len(res.GetLeaves()) == 2)
401 ks = res.GetLeaves().keys()
402 self.assertTrue('[*]n1cccc1' in ks)
403 self.assertTrue('[*]NC1CO1' in ks)
404
405 m = Chem.MolFromSmiles('c1cccn1C(=O)n1c(C)ccc1')
406 res = RecapDecompose(m, onlyUseReactions=[0])
407 self.assertTrue(res)
408 self.assertTrue(len(res.GetLeaves()) == 2)
409 ks = res.GetLeaves().keys()
410 self.assertTrue('[*]n1cccc1C' in ks)
411
413 m = Chem.MolFromSmiles('C1CC1N(C1NC1)C1OC1')
414 res = RecapDecompose(m)
415 self.assertTrue(res)
416 self.assertTrue(len(res.GetLeaves()) == 3)
417 ks = res.GetLeaves().keys()
418 self.assertTrue('[*]C1CC1' in ks)
419 self.assertTrue('[*]C1CO1' in ks)
420 self.assertTrue('[*]C1CN1' in ks)
421
422 m = Chem.MolFromSmiles('c1ccccc1N(C1NC1)C1OC1')
423 res = RecapDecompose(m)
424 self.assertTrue(res)
425 self.assertTrue(len(res.GetLeaves()) == 3)
426 ks = res.GetLeaves().keys()
427 self.assertTrue('[*]c1ccccc1' in ks)
428 self.assertTrue('[*]C1CO1' in ks)
429 self.assertTrue('[*]C1CN1' in ks)
430
431 m = Chem.MolFromSmiles('c1ccccc1N(c1ncccc1)C1OC1')
432 res = RecapDecompose(m)
433 self.assertTrue(res)
434 self.assertTrue(len(res.GetLeaves()) == 3)
435 ks = res.GetLeaves().keys()
436 self.assertTrue('[*]c1ccccc1' in ks)
437 self.assertTrue('[*]c1ccccn1' in ks)
438 self.assertTrue('[*]C1CO1' in ks)
439
440 m = Chem.MolFromSmiles('c1ccccc1N(c1ncccc1)c1ccco1')
441 res = RecapDecompose(m)
442 self.assertTrue(res)
443 self.assertTrue(len(res.GetLeaves()) == 3)
444 ks = res.GetLeaves().keys()
445 self.assertTrue('[*]c1ccccc1' in ks)
446 self.assertTrue('[*]c1ccccn1' in ks)
447 self.assertTrue('[*]c1ccco1' in ks)
448
449 m = Chem.MolFromSmiles('C1CCCCN1C1CC1')
450 res = RecapDecompose(m)
451 self.assertTrue(res)
452 self.assertTrue(len(res.GetLeaves()) == 2)
453 ks = res.GetLeaves().keys()
454 self.assertTrue('[*]N1CCCCC1' in ks)
455 self.assertTrue('[*]C1CC1' in ks)
456
457 m = Chem.MolFromSmiles('C1CCC2N1CC2')
458 res = RecapDecompose(m)
459 self.assertTrue(res)
460 self.assertTrue(len(res.GetLeaves()) == 0)
461
463 m = Chem.MolFromSmiles('C1CC1OC1OC1')
464 res = RecapDecompose(m)
465 self.assertTrue(res)
466 self.assertTrue(len(res.GetLeaves()) == 2)
467 ks = res.GetLeaves().keys()
468 self.assertTrue('[*]C1CC1' in ks)
469 self.assertTrue('[*]C1CO1' in ks)
470
471 m = Chem.MolFromSmiles('C1CCCCO1')
472 res = RecapDecompose(m)
473 self.assertTrue(res)
474 self.assertTrue(len(res.GetLeaves()) == 0)
475
476 m = Chem.MolFromSmiles('c1ccccc1OC1OC1')
477 res = RecapDecompose(m)
478 self.assertTrue(res)
479 self.assertTrue(len(res.GetLeaves()) == 2)
480 ks = res.GetLeaves().keys()
481 self.assertTrue('[*]c1ccccc1' in ks)
482 self.assertTrue('[*]C1CO1' in ks)
483
484 m = Chem.MolFromSmiles('c1ccccc1Oc1ncccc1')
485 res = RecapDecompose(m)
486 self.assertTrue(res)
487 self.assertTrue(len(res.GetLeaves()) == 2)
488 ks = res.GetLeaves().keys()
489 self.assertTrue('[*]c1ccccc1' in ks)
490 self.assertTrue('[*]c1ccccn1' in ks)
491
493 m = Chem.MolFromSmiles('ClC=CBr')
494 res = RecapDecompose(m)
495 self.assertTrue(res)
496 self.assertTrue(len(res.GetLeaves()) == 2)
497 ks = res.GetLeaves().keys()
498 self.assertTrue('[*]CCl' in ks)
499 self.assertTrue('[*]CBr' in ks)
500
501 m = Chem.MolFromSmiles('C1CC=CC1')
502 res = RecapDecompose(m)
503 self.assertTrue(res)
504 self.assertTrue(len(res.GetLeaves()) == 0)
505
507 m = Chem.MolFromSmiles('c1cccn1CCCC')
508 res = RecapDecompose(m)
509 self.assertTrue(res)
510 self.assertTrue(len(res.GetLeaves()) == 2)
511 ks = res.GetLeaves().keys()
512 self.assertTrue('[*]n1cccc1' in ks)
513 self.assertTrue('[*]CCCC' in ks)
514
515 m = Chem.MolFromSmiles('c1ccc2n1CCCC2')
516 res = RecapDecompose(m)
517 self.assertTrue(res)
518 self.assertTrue(len(res.GetLeaves()) == 0)
519
521 m = Chem.MolFromSmiles('C1CC(=O)N1CCCC')
522 res = RecapDecompose(m, onlyUseReactions=[8])
523 self.assertTrue(res)
524 self.assertTrue(len(res.GetLeaves()) == 2)
525 ks = res.GetLeaves().keys()
526 self.assertTrue('[*]N1CCC1=O' in ks)
527 self.assertTrue('[*]CCCC' in ks)
528
529 m = Chem.MolFromSmiles('O=C1CC2N1CCCC2')
530 res = RecapDecompose(m)
531 self.assertTrue(res)
532 self.assertTrue(len(res.GetLeaves()) == 0)
533
535 m = Chem.MolFromSmiles('c1ccccc1c1ncccc1')
536 res = RecapDecompose(m)
537 self.assertTrue(res)
538 self.assertTrue(len(res.GetLeaves()) == 2)
539 ks = res.GetLeaves().keys()
540 self.assertTrue('[*]c1ccccc1' in ks)
541 self.assertTrue('[*]c1ccccn1' in ks)
542
543 m = Chem.MolFromSmiles('c1ccccc1C1CC1')
544 res = RecapDecompose(m)
545 self.assertTrue(res)
546 self.assertTrue(len(res.GetLeaves()) == 0)
547
549 m = Chem.MolFromSmiles('c1cccn1c1ccccc1')
550 res = RecapDecompose(m)
551 self.assertTrue(res)
552 self.assertTrue(len(res.GetLeaves()) == 2)
553 ks = res.GetLeaves().keys()
554 self.assertTrue('[*]n1cccc1' in ks)
555 self.assertTrue('[*]c1ccccc1' in ks)
556
558 m = Chem.MolFromSmiles('CCCNS(=O)(=O)CC')
559 res = RecapDecompose(m)
560 self.assertTrue(res)
561 self.assertTrue(len(res.GetLeaves()) == 2)
562 ks = res.GetLeaves().keys()
563 self.assertTrue('[*]NCCC' in ks)
564 self.assertTrue('[*]S(=O)(=O)CC' in ks)
565
566 m = Chem.MolFromSmiles('c1cccn1S(=O)(=O)CC')
567 res = RecapDecompose(m)
568 self.assertTrue(res)
569 self.assertTrue(len(res.GetLeaves()) == 2)
570 ks = res.GetLeaves().keys()
571 self.assertTrue('[*]n1cccc1' in ks)
572 self.assertTrue('[*]S(=O)(=O)CC' in ks)
573
574 m = Chem.MolFromSmiles('C1CNS(=O)(=O)CC1')
575 res = RecapDecompose(m)
576 self.assertTrue(res)
577 self.assertTrue(len(res.GetLeaves()) == 0)
578
580 m = Chem.MolFromSmiles('c1ccccc1n1cccc1')
581 res = RecapDecompose(m)
582 self.assertTrue(res)
583 self.assertTrue(len(res.GetLeaves()) == 2)
584 m = Chem.MolFromSmiles('c1ccccc1[n+]1ccccc1')
585 res = RecapDecompose(m)
586 self.assertTrue(res)
587 self.assertTrue(len(res.GetLeaves()) == 0)
588
589 m = Chem.MolFromSmiles('C1CC1NC(=O)CC')
590 res = RecapDecompose(m)
591 self.assertTrue(res)
592 self.assertTrue(len(res.GetLeaves()) == 2)
593 m = Chem.MolFromSmiles('C1CC1[NH+]C(=O)CC')
594 res = RecapDecompose(m)
595 self.assertTrue(res)
596 self.assertTrue(len(res.GetLeaves()) == 0)
597
598 m = Chem.MolFromSmiles('C1CC1NC(=O)NC1CCC1')
599 res = RecapDecompose(m)
600 self.assertTrue(res)
601 self.assertTrue(len(res.GetLeaves()) == 2)
602 m = Chem.MolFromSmiles('C1CC1[NH+]C(=O)[NH+]C1CCC1')
603 res = RecapDecompose(m)
604 self.assertTrue(res)
605 self.assertTrue(len(res.GetLeaves()) == 0)
606
607 unittest.main()
608