Mercurial > hg > python
annotate cluster.py @ 76:cba9cb3a6797
from python3.13 with block boundary tracking added by HST
| author | Henry S Thompson <ht@inf.ed.ac.uk> |
|---|---|
| date | Thu, 12 Mar 2026 18:06:10 +0000 |
| parents | 668d788cac47 |
| children |
| rev | line source |
|---|---|
| 23 | 1 #!/usr/bin/python3 |
| 2 '''Read a minutor block lockation tsv and resort it to show clusters''' | |
| 3 import sys | |
| 4 | |
| 5 from util import * | |
| 6 | |
| 7 if len(sys.argv)==1 or sys.argv[1]=='-h': | |
| 8 print("""Usage: cluster.py [-h] [-c n] [-s] infile.tsv [outfile.tsv] | |
| 9 n is cluster diameter, default is 5 | |
| 10 -s for strict circle distance, reduces number of candidates | |
| 11 default outfile is [infile]_c[n].tsv""") | |
| 12 exit(1) | |
| 13 if sys.argv[1]=='-c': | |
| 14 sys.argv.pop(1) | |
| 15 n=float(sys.argv.pop(1)) | |
| 16 else: | |
| 17 n=5.0 | |
| 18 | |
| 19 if sys.argv[1]=='-s': | |
| 20 sys.argv.pop(1) | |
| 21 strict=True | |
| 22 else: | |
| 23 strict=False | |
| 24 | |
| 25 infile_name=sys.argv.pop(1) | |
| 26 if len(sys.argv)>1: | |
| 27 outfile_name=sys.argv.pop(1) | |
| 28 else: | |
| 29 outfile_name="%s_c%s.tsv"%(infile_name.split('.')[0],n) | |
| 30 | |
| 31 cc=[] | |
| 32 | |
| 33 with open(infile_name,'r') as infile: | |
| 34 with open(outfile_name,'w') as outfile: | |
| 35 l=infile.readline().rstrip() | |
| 36 print(l,file=outfile) | |
| 37 ff=PPAT.split(l) | |
| 38 (nr,ox,oy,oz)=intsMaybe(ff) | |
| 39 home=(float(ox),float(oy),float(oz)) | |
| 40 et=ff[9] | |
| 41 l=infile.readline().rstrip() | |
| 42 print(l,file=outfile) | |
|
34
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
43 l2=PPAT.split(l) |
|
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
44 if l2[-1]=='bounded': |
|
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
45 orad=int(l2[1]) |
|
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
46 ymin=ymax=-1 |
|
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
47 else: |
|
668d788cac47
handle unbounded Y, landing below target elevation
Henry Thompson <ht@markup.co.uk>
parents:
23
diff
changeset
|
48 (orad,ymin,ymax)=intsMaybe(l2) |
| 23 | 49 _=infile.readline() |
| 50 for l in infile: | |
| 51 found=False | |
| 52 (_,td,qd)=l.rstrip().split('\t') | |
| 53 q=[float(i) for i in qd.split(',')] | |
| 54 td=float(td) | |
| 55 if strict and td>orad: | |
| 56 # Enforce a circular region, not square | |
| 57 print(td,orad,file=sys.stderr) | |
| 58 break | |
| 59 for c in cc: | |
| 60 for p in c: | |
| 61 if d(p,q)<=n: | |
| 62 c.append(q) | |
| 63 found=True | |
| 64 break | |
| 65 if found: | |
| 66 break | |
| 67 if not found: | |
| 68 cc.append([q]) | |
| 69 oc=cc | |
| 70 cc=[] # lose | |
| 71 w=0 | |
| 72 ow=-1 | |
| 73 nc=[] # win | |
| 74 while True: | |
| 75 for i,c in enumerate(oc): | |
| 76 win=False | |
| 77 for p in c: | |
| 78 for g in oc[i+1:]: | |
| 79 for q in g: | |
| 80 if d(p,q)<=n: | |
| 81 win=True | |
| 82 w+=1 | |
| 83 nc.append(c+g) | |
| 84 break | |
| 85 if win: | |
| 86 break | |
| 87 if win: | |
| 88 break | |
| 89 if not win: | |
| 90 cc.append(c) | |
| 91 print(len(cc),len(nc),ow,w,file=sys.stderr) | |
| 92 if ow==w: | |
| 93 break | |
| 94 ow=w | |
| 95 oc=nc | |
| 96 nc=[] | |
| 97 for c in sorted(cc,reverse=True,key=lambda x:len(x)): | |
| 98 print(c,file=outfile) | |
| 99 | |
| 100 | |
| 101 | |
| 102 |