-
Notifications
You must be signed in to change notification settings - Fork 5
Expand file tree
/
Copy pathfixstardirs.py
More file actions
981 lines (885 loc) · 37.5 KB
/
fixstardirs.py
File metadata and controls
981 lines (885 loc) · 37.5 KB
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
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
# -*- coding: utf-8 -*-
import math
import datetime
import sweastrology as swe
import os
import houses
import common # Morinus의 ephe 경로 사용 (common.common.ephepath)
import chart
import fixstars
import primdirs
import astrology
import common # Morinus의 ephe 경로 사용 (common.common.ephepath)
#try: swe.set_ephe_path(common.common.ephepath)
#except: pass
def _ensure_swisseph_path():
# 1) Morinus 설정에서 우선 시도
cands = []
try:
ep = getattr(common.common, 'ephepath', '')
if ep:
ep = os.path.expandvars(os.path.expanduser(ep))
if not os.path.isabs(ep):
# 모듈 파일 위치 기준으로 절대경로 변환
base = os.path.abspath(os.path.dirname(__file__))
ep_abs = os.path.normpath(os.path.join(base, ep))
cands += [ep_abs]
cands += [ep, os.path.join(ep, 'SWEP', 'Ephem')]
except Exception:
pass
# 2) 모듈 기준 상대 경로들
here = os.path.abspath(os.path.dirname(__file__))
cands += [
os.path.join(here, 'SWEP', 'Ephem'),
os.path.join(here, 'Ephem'),
here,
os.environ.get('SWEPH_PATH', '') # 환경변수 지원
]
# 3) 실제 Swiss Ephemeris 바이너리(se*.se1)가 있는 첫 경로로 설정
for p in cands:
if not p:
continue
p = os.path.normpath(p)
try:
if os.path.isdir(p) and any(fn.lower().startswith(('sepl_', 'sef', 'semo', 'sea'))
for fn in os.listdir(p)):
swe.swe_set_ephe_path(p)
return
except Exception:
continue
# 여기까지 못 잡았으면 그대로 두되, 이후 호출부에서 다시 시도/에러 노출
# (raise로 끊고 싶으면 RuntimeError(...)를 던지도록 바꿔도 됨)
# 최초 1회 시도(실패해도 지나가게)
try:
_ensure_swisseph_path()
except Exception:
pass
import re
NAIBOD_COEFF = primdirs.PrimDirs.staticData[primdirs.PrimDirs.NAIBOD][3] # 1.01456164
# [ADD] ---- Primary Directions Key를 따르는 arc→years 변환 헬퍼 ----
DEG = math.pi / 180.0
TROPICAL_YEAR = 365.2421897
def _norm360(x):
return x % 360.0
def _sun_coord_deg(jd_ut, equatorial=False):
_ensure_swisseph_path()
"""equatorial=True면 적경(°), False면 황경(°) 반환"""
flags = astrology.SEFLG_SWIEPH | (astrology.SEFLG_EQUATORIAL if equatorial else 0)
try:
vals, _ = swe.swe_calc_ut(jd_ut, astrology.SE_SUN, flags)
except Exception:
# 파일(ephe) 없으면 내부 알고리즘(MOSEPH)로 재시도 → 파일 불필요
flags = astrology.SEFLG_MOSEPH | (astrology.SEFLG_EQUATORIAL if equatorial else 0)
vals, _ = swe.swe_calc_ut(jd_ut, astrology.SE_SUN, flags)
lon = vals[0]
return _norm360(lon)
def _true_solar_arc_years(horoscope, options, arc_deg, direct):
"""진태양(적경/황경) 동적키: 출생 시점에서 태양 좌표가 arc_deg만큼 변하는 시점까지의 '연수' 반환"""
equ = (options.pdkeyd == primdirs.PrimDirs.TRUESOLAREQUATORIALARC)
jd0 = _birth_jd_ut(horoscope)
base = _sun_coord_deg(jd0, equ)
# Direct이거나(+) / Converse에서도 'Use regressive'가 꺼져 있으면(+) 정방으로, 그렇지 않으면 역방(-)
sign = +1.0 if (direct or not getattr(options, 'useregressive', False)) else -1.0
# 하루 이동량으로 초기 추정: 적경/황경 각각에 대해 실제 하루차로 계측
if sign > 0:
step = _norm360(_sun_coord_deg(jd0 + 1.0, equ) - base)
else:
step = _norm360(base - _sun_coord_deg(jd0 - 1.0, equ))
if step < 1e-6:
step = 0.985647 # 안전한 기본치(황경 기준 하루 이동량)
# wrap-around로 2° 근처가 잡히는 드문 케이스 보정(태양은 ~1°/day가 정상)
if step > 2.0:
step = 360.0 - step
t0 = 0.0
t1 = sign * (arc_deg / step) # days
MAX_DAYS = 400.0 # 1년(≈365d)보다 약간 큰 안전 범위
if abs(t1) > MAX_DAYS:
t1 = math.copysign(MAX_DAYS, t1)
def f(t):
pos = _sun_coord_deg(jd0 + t, equ)
if sign > 0:
return _norm360(pos - base) - arc_deg
else:
return _norm360(base - pos) - arc_deg
f0 = f(t0)
f1 = f(t1)
for _ in range(8):
denom = (f1 - f0)
if abs(denom) < 1e-12:
# 기울기 소실 → 폭주 방지: 이분법성 수렴으로 전환
t2 = 0.5 * (t1 + t0)
else:
t2 = t1 - f1 * (t1 - t0) / denom
# 비정상 추정치(무한/NaN/범위초과) 클램프
if not math.isfinite(t2) or abs(t2) > MAX_DAYS:
t2 = 0.5 * (t1 + t0)
t0, f0, t1, f1 = t1, f1, t2, f(t2)
if abs(t1 - t0) < 1e-6 or abs(f1) < 1e-9:
break
return abs(t1)
def _birthday_solar_arc_years(horoscope, options, arc_deg):
"""생일태양(적경/황경) 동적키: 출생일과 다음날 태양 이동량으로 환산"""
equ = (options.pdkeyd == primdirs.PrimDirs.BIRTHDAYSOLAREQUATORIALARC)
jd0 = _birth_jd_ut(horoscope)
c0 = _sun_coord_deg(jd0, equ)
c1 = _sun_coord_deg(jd0 + 1.0, equ)
d = _norm360(c1 - c0)
return (arc_deg / d) if d != 0.0 else 0.0
def _arc_to_years_from_primary_key(horoscope, options, arc_deg, direct):
"""PrimDirs의 calcTime 로직을 최소 구현(정적키/동적키 모두 지원)"""
if options.pdkeydyn:
# 동적 키
if options.pdkeyd in (primdirs.PrimDirs.TRUESOLAREQUATORIALARC,
primdirs.PrimDirs.TRUESOLARECLIPTICALARC):
return _true_solar_arc_years(horoscope, options, arc_deg, direct)
else:
return _birthday_solar_arc_years(horoscope, options, arc_deg)
else:
# 정적 키
if options.pdkeys == primdirs.PrimDirs.CUSTOMER:
val = (options.pdkeydeg + options.pdkeymin/60.0 + options.pdkeysec/3600.0) # deg/year
return (arc_deg / val) if val != 0.0 else 0.0
else:
coeff = primdirs.PrimDirs.staticData[options.pdkeys][primdirs.PrimDirs.COEFF] # years/deg
return arc_deg * coeff
# [END ADD]
_FIXSTAR_CAT_DB = None
def _adlat(phi_deg, dec_deg):
"""ADlat = asin(tan φ * tan δ). |tanφ·tanδ| > 1 → ASC/DSC 부재."""
val = math.tan(math.radians(phi_deg)) * math.tan(math.radians(dec_deg))
if abs(val) > 1.0:
return None
return math.degrees(math.asin(val))
def _ramc_pack(horoscope):
"""RAMC, RAIC, AOASC, DODESC (모리누스와 동일 정의)"""
try:
ramc = horoscope.houses.ascmc2[houses.Houses.MC][houses.Houses.RA]
except Exception:
ramc = _get_ramc0_deg(horoscope)
raic = (ramc + 180.0) % 360.0
aoasc = (ramc + 90.0) % 360.0
dodesc = (ramc + 270.0) % 360.0
return ramc, raic, aoasc, dodesc
def _arc_direct(base_deg, target_deg):
"""직행(D): base→target 으로 전진해야 하는 양(0..360)"""
return (target_deg - base_deg) % 360.0
def _arc_converse(base_deg, target_deg):
"""역행(C): base←target 으로 후진해야 하는 양(0..360)"""
return (base_deg - target_deg) % 360.0
def _fixstars_cat_paths():
# 네가 말한 경로: 작업폴더\SWEP\Ephem 도 같이 스캔
import os
paths = []
try:
ep = getattr(common.common, 'ephepath', '')
if ep:
paths += [
os.path.join(ep, 'sefstars.txt'),
os.path.join(ep, 'SWEP', 'Ephem', 'sefstars.txt'),
os.path.join(ep, 'fixstars.cat'),
os.path.join(ep, 'fixedstars.cat'),
os.path.join(ep, 'SWEP', 'Ephem', 'fixstars.cat'),
]
except Exception:
pass
here = os.path.dirname(__file__)
paths += [
os.path.join(here, 'SWEP', 'Ephem', 'sefstars.txt'),
os.path.join(here, 'sefstars.txt'),
os.path.join(here, 'SWEP', 'Ephem', 'fixstars.cat'),
os.path.join(here, 'fixstars.cat'),
os.path.join(here, 'fixedstars.cat'),
]
# 존재하는 것만
return [p for p in paths if os.path.isfile(p)]
def _load_fixstars_cat():
"""
fixstars/fixedstars.cat을 1회 캐시 로드:
code -> {name, ra_j2000_deg, dec_j2000_deg, pm_ra_sec, pm_dec_arcsec}
(콤마 구분 포맷 기준. 없으면 0으로 간주)
"""
global _FIXSTAR_CAT_DB
if _FIXSTAR_CAT_DB is not None:
return _FIXSTAR_CAT_DB
db = {}
for path in _fixstars_cat_paths():
try:
f = open(path, 'r')
for line in f:
s = line.strip()
if not s or s.startswith('#'):
continue
parts = [p.strip() for p in s.split(',')]
if len(parts) < 9:
continue
name = parts[0]
code = parts[1].lstrip(',')
# RA
try:
ra_h = float(parts[3]); ra_m = float(parts[4]); ra_s = float(parts[5])
ra_deg = (ra_h + ra_m/60.0 + ra_s/3600.0) * 15.0
except:
continue
# Dec
try:
dec_d_raw = parts[6]
dec_d = abs(float(dec_d_raw))
dec_m = float(parts[7]); dec_s = float(parts[8])
sign = -1.0 if str(dec_d_raw).strip().startswith('-') else 1.0
dec_deg = sign * (dec_d + dec_m/60.0 + dec_s/3600.0)
except:
continue
pm_ra_sec = 0.0
pm_dec_arcsec = 0.0
try:
# 파일별 단위 판별: sefstars.txt( mas/yr, RA에 cosδ 포함 ) vs legacy( 초/세기, 호초/세기 )
is_sef = os.path.basename(path).lower().startswith('sefstars')
if is_sef:
# parts[9]=pmRA(mas/yr * cosδ0), parts[10]=pmDE(mas/yr)
pm_ra_masy = float(parts[9]) if len(parts) >= 10 else 0.0
pm_de_masy = float(parts[10]) if len(parts) >= 11 else 0.0
# 내부 형식으로 변환:
# RA: (mas/yr * cosδ0) → 초(시간)/세기 = ( (mas/yr)/cosδ0 ) / 150
# Dec: mas/yr → 호초/세기 = (mas/yr) / 10
cosd = max(1e-12, math.cos(math.radians(dec_deg)))
pm_ra_sec = (pm_ra_masy / cosd) / 150.0
pm_dec_arcsec = pm_de_masy / 10.0
elif len(parts) >= 11:
pm_ra_sec = float(parts[9])
pm_dec_arcsec = float(parts[10])
except:
pm_ra_sec = 0.0; pm_dec_arcsec = 0.0
db[code] = {
'name': name,
'ra_j2000_deg': ra_deg,
'dec_j2000_deg': dec_deg,
'pm_ra_sec': pm_ra_sec,
'pm_dec_arcsec': pm_dec_arcsec
}
f.close()
except Exception:
continue
_FIXSTAR_CAT_DB = db
return _FIXSTAR_CAT_DB
# --- Python 2/3 unicode helpers ---
try:
unicode # Py2
except NameError:
unicode = str
basestring = (str, )
# --- fixed stars catalog path helper (supports both names) ---
def _fixstars_cat_path():
base = getattr(common.common, 'ephepath', '')
p0 = os.path.join(base, 'sefstars.txt')
p1 = os.path.join(base, 'SWEP', 'Ephem', 'sefstars.txt')
p2 = os.path.join(base, 'fixstars.cat') # legacy
p3 = os.path.join(base, 'SWEP', 'Ephem', 'fixstars.cat') # legacy in SWEP\Ephem
p4 = os.path.join(base, 'fixedstars.cat') # legacy alt
if os.path.exists(p0): return p0 # 1) ep\sefstars.txt
if os.path.exists(p1): return p1 # 2) ep\SWEP\Ephem\sefstars.txt
if os.path.exists(p2): return p2 # 3) legacy
if os.path.exists(p3): return p3 # 4) legacy in SWEP\Ephem
if os.path.exists(p4): return p4 # 5) legacy alt
return None
def _cat_all_names_in_order():
"""
카탈로그의 '첫 컬럼 이름'을 위에서 아래로 그대로 읽어 반환.
콤마/세미콜론/공백 포맷 모두 지원.
"""
path = _fixstars_cat_path()
if not path:
return []
names = []
f = open(path, "r")
try:
for ln in f:
if not isinstance(ln, unicode):
try: ln = ln.decode('utf-8', 'ignore')
except: pass
s = ln.strip()
if (not s) or s.startswith("#"):
continue
if "," in s:
nm = s.split(",", 1)[0].strip()
elif ";" in s:
nm = s.split(";", 1)[0].strip()
else:
nm = s.split()[0].strip()
if nm:
names.append(nm if isinstance(nm, unicode) else unicode(nm))
finally:
try: f.close()
except: pass
return names
def _ra_dec_star_ofdate_from_code(code, jd_ut):
"""
code(예: 'alTau')로 카탈로그 DB에서 J2000 좌표+고유운동을 얻고,
of-date로 세차해서 RA/Dec(deg)를 반환.
"""
db = _load_fixstars_cat()
if code not in db:
raise RuntimeError(u"fixstars.cat에 항성 코드가 없습니다: %s" % code)
ra0 = db[code]['ra_j2000_deg']
de0 = db[code]['dec_j2000_deg']
pmra = db[code].get('pm_ra_sec', 0.0)
pmde = db[code].get('pm_dec_arcsec', 0.0)
# (이미 네 파일에 있는) 고유운동/세차 유틸 사용
# 고유운동: 세기→연 환산 반영
ra1, de1 = _apply_proper_motion_j2000(ra0, de0, jd_ut, pmra, pmde)
# TT = UT + ΔT
jd_tt = jd_ut + swe.swe_deltat(jd_ut)
ra2, de2 = _j2000_to_ofdate(ra1, de1, jd_tt)
return ra2, de2, db[code].get('name', code)
# ===== B1950(FK4) → J2000(FK5) → of-date(출생시) 변환 유틸 =====
J2000_JD = 2451545.0 # TT
def _apply_proper_motion_j2000(ra_deg, dec_deg, jd_ut, pm_ra_sec_century, pm_dec_arcsec_century):
"""
ra/dec: J2000 기준(도)
pm_ra_sec_century: RA 고유운동 [초(시간)/세기]
pm_dec_arcsec_century: Dec 고유운동 [호초/세기]
jd_ut: 출생 UTJulianDay
"""
# J2000 으로부터 경과 '연'
dt_years = (jd_ut - J2000_JD) / 365.2421897
# *** 세기 → 연 으로 환산 ***
# RA: 초(시간)/세기 → (초/연) → 도/연
pm_ra_sec_year = float(pm_ra_sec_century) / 100.0
d_ra_deg_per_yr = pm_ra_sec_year * (15.0 / 3600.0) # 1s(time)=15″
# Dec: 호초/세기 → (호초/연) → 도/연
pm_dec_arcsec_year = float(pm_dec_arcsec_century) / 100.0
d_dec_deg_per_yr = pm_dec_arcsec_year / 3600.0
ra_new = (ra_deg + d_ra_deg_per_yr * dt_years) % 360.0
dec_new = dec_deg + d_dec_deg_per_yr * dt_years
if dec_new > 90.0: dec_new = 90.0
elif dec_new < -90.0: dec_new = -90.0
return ra_new, dec_new
def _rad(d): return d * (math.pi/180.0)
def _deg(r): return r * (180.0/math.pi)
def _ra_dec_to_vec(ra_deg, dec_deg):
ra = _rad(ra_deg); dec = _rad(dec_deg)
cosd = math.cos(dec)
return [math.cos(ra)*cosd, math.sin(ra)*cosd, math.sin(dec)]
def _vec_to_ra_dec(v):
x,y,z = v
rxy = math.hypot(x,y)
ra = math.atan2(y, x)
if ra < 0.0: ra += 2.0*math.pi
dec = math.atan2(z, rxy)
return _deg(ra), _deg(dec)
def _norm_vec(v):
x,y,z = v
n = math.sqrt(x*x+y*y+z*z)
if n == 0.0: return [0.0,0.0,0.0]
return [x/n, y/n, z/n]
def _mat_vec(M, v):
return [
M[0][0]*v[0] + M[0][1]*v[1] + M[0][2]*v[2],
M[1][0]*v[0] + M[1][1]*v[1] + M[1][2]*v[2],
M[2][0]*v[0] + M[2][1]*v[1] + M[2][2]*v[2],
]
# FK4(B1950) → FK5(J2000) (Aoki et al. 1983 근사; proper motion 미사용)
_M_FK4_TO_FK5 = [
[0.9999256782, -0.0111820611, -0.0048579477],
[0.0111820610, 0.9999374784, -0.0000271765],
[0.0048579479, -0.0000271474, 0.9999881997],
]
_A_FK4_TO_FK5 = [-1.62557e-6, -0.31919e-6, -0.13843e-6] # E-terms 보정 상수
def _fk4_b1950_to_fk5_j2000(ra50_deg, dec50_deg):
r = _ra_dec_to_vec(ra50_deg, dec50_deg)
# r_J2000 ≈ M*(r_FK4) + A (정규화)
rj = _mat_vec(_M_FK4_TO_FK5, r)
rj = [rj[0] + _A_FK4_TO_FK5[0],
rj[1] + _A_FK4_TO_FK5[1],
rj[2] + _A_FK4_TO_FK5[2]]
rj = _norm_vec(rj)
return _vec_to_ra_dec(rj)
# IAU 1976 (Lieske 1979) 프리세션: J2000 → of-date (TT)
def _precession_matrix_J2000_to_date(jd_tt):
t = (jd_tt - J2000_JD) / 36525.0 # Julian centuries from J2000
# 호(arcsec)
zeta = (2306.2181*t + 0.30188*t*t + 0.017998*t*t*t)
z = (2306.2181*t + 1.09468*t*t + 0.018203*t*t*t)
theta = (2004.3109*t - 0.42665*t*t - 0.041833*t*t*t)
zeta = _rad(zeta/3600.0); z = _rad(z/3600.0); theta = _rad(theta/3600.0)
cz, sz = math.cos(z), math.sin(z)
czeta, szeta = math.cos(zeta), math.sin(zeta)
cth, sth = math.cos(theta), math.sin(theta)
# Rz(-z) * Ry(theta) * Rz(-zeta)
Rz1 = [[ cz, sz, 0], [-sz, cz, 0], [0,0,1]]
Ry = [[ cth, 0, -sth], [0,1,0], [ sth,0, cth]]
Rz2 = [[ czeta, szeta, 0], [-szeta, czeta, 0], [0,0,1]]
# multiply: R = Rz1 * Ry * Rz2
def _mm(A,B):
return [[A[0][0]*B[0][0]+A[0][1]*B[1][0]+A[0][2]*B[2][0],
A[0][0]*B[0][1]+A[0][1]*B[1][1]+A[0][2]*B[2][1],
A[0][0]*B[0][2]+A[0][1]*B[1][2]+A[0][2]*B[2][2]],
[A[1][0]*B[0][0]+A[1][1]*B[1][0]+A[1][2]*B[2][0],
A[1][0]*B[0][1]+A[1][1]*B[1][1]+A[1][2]*B[2][1],
A[1][0]*B[0][2]+A[1][1]*B[1][2]+A[1][2]*B[2][2]],
[A[2][0]*B[0][0]+A[2][1]*B[1][0]+A[2][2]*B[2][0],
A[2][0]*B[0][1]+A[2][1]*B[1][1]+A[2][2]*B[2][1],
A[2][0]*B[0][2]+A[2][1]*B[1][2]+A[2][2]*B[2][2]]]
return _mm(_mm(Rz1, Ry), Rz2)
def _j2000_to_ofdate(ra2000_deg, dec2000_deg, jd_tt):
R = _precession_matrix_J2000_to_date(jd_tt)
v0 = _ra_dec_to_vec(ra2000_deg, dec2000_deg)
v = _mat_vec(R, v0)
v = _norm_vec(v)
return _vec_to_ra_dec(v)
def _clamp(x, a, b):
return a if x < a else b if x > b else x
def _norm360(x):
x = x % 360.0
return x if x >= 0 else x + 360.0
def asc_dsc_exists(phi_deg, dec_deg):
# |tan φ * tan δ| <= 1 이어야 상승/하강 존재
x = abs(math.tan(phi_deg*DEG) * math.tan(dec_deg*DEG))
return x <= 1.0 - 1e-12 # 임계 근접 수치 폭주 방지 마진
def _pair_arc_and_label(diff_deg):
"""
primdirs.create()와 같은 규칙:
- 작은호(≤180)를 기본으로 취하고, 그때의 레이블이 D(직행)
- 보완호(360-작은호)는 C(역행)
- diff가 180초과면 레이블을 뒤집는 효과가 발생
"""
d = diff_deg % 360.0
if d <= 180.0:
return (d, 'D'), (360.0 - d, 'C')
else:
# 작은호는 360-d (역행이 더 가깝다)
return (360.0 - d, 'C'), (d, 'D')
def primary_arcs(ramc0_deg, ramc_evt_deg):
# Direct/Converse 모두 "미래 나이"로 매핑
arc_dir = _norm360(ramc_evt_deg - ramc0_deg)
arc_conv = _norm360(ramc0_deg - ramc_evt_deg)
return arc_dir, arc_conv
def arc_to_age_years_naibod(arc_deg):
return arc_deg / 360.0
def _get_ramc0_deg(horoscope):
# 1) 이미 계산된 RAMC 보유 시 재사용
for cand in ("ramc", "RAMC", "ramc_deg", "ramc0"):
if hasattr(horoscope, cand):
val = getattr(horoscope, cand)
try:
return float(val) % 360.0
except:
pass
# 2) Fallback: LST(그리니치 항성시 + 경도)로 계산
jd_ut = float(getattr(horoscope, "jd_ut", getattr(horoscope, "jdut", 0.0)))
lon_deg = float(getattr(horoscope, "lon", getattr(horoscope, "longitude", 0.0)))
sid_hours = swe.swe_sidtime(jd_ut) # Greenwich sidereal time [hours]
lst_deg = (sid_hours * 15.0 + lon_deg) % 360.0
return lst_deg
def _birth_jd_ut(horoscope):
try:
return float(horoscope.time.jd)
except:
return float(getattr(horoscope, "jd_ut", getattr(horoscope, "jdut", 0.0)))
def _observer_lat(horoscope):
try:
return float(horoscope.place.lat)
except:
pass
try:
deglat = float(horoscope.place.deglat)
minlat = float(getattr(horoscope.place, 'minlat', 0.0))
seclat = float(getattr(horoscope.place, 'seclat', 0.0))
north = bool(getattr(horoscope.place, 'north', True))
lat = abs(deglat) + minlat/60.0 + seclat/3600.0
if not north:
lat = -lat
return lat
except:
return float(getattr(horoscope, "lat", getattr(horoscope, "latitude", 0.0)))
def _calendar_flag(chrt, options):
# 1=Gregorian, 0=Julian
try:
calobj = getattr(chrt, "time", None)
if calobj is not None:
return 0 if calobj.cal == chart.Time.JULIAN else 1
except Exception:
pass
cal = getattr(options, "calendar", "greg")
if isinstance(cal, int):
return 1 if cal != 0 else 0
cal = str(cal).lower()
return 1 if ("greg" in cal or cal == "g") else 0
def _cat_name_to_code_map():
"""카탈로그 DB에서 name→code 매핑 생성."""
db = _load_fixstars_cat()
m = {}
for code, rec in db.items():
nm = rec.get('name', u'')
if not isinstance(nm, unicode):
try: nm = unicode(nm)
except: pass
m[nm] = code
return m
def _selected_star_codes(options):
"""
우선순위:
1) options.fixstars 가 dict이면 그 'key(코드)'만 사용
2) options.pdfixstarssel 이 dict이면 True인 '이름'을 code로 매핑해 사용
(둘 다 없으면 빈 리스트)
반환: ['alTau', 'spica', ...] 처럼 카탈로그 'code' 리스트
"""
# 1) fixstars: {code: ...}
fs = getattr(options, 'fixstars', None)
if isinstance(fs, dict) and len(fs) > 0:
try:
keys_iter = fs.iterkeys()
except AttributeError:
keys_iter = fs.keys()
out = []
for k in keys_iter:
ku = k if isinstance(k, unicode) else unicode(k)
ku = ku.strip().lstrip(',') # 혹시 앞에 콤마가 붙어 있다면 제거
if ku:
out.append(ku)
if out:
return out
# 2) pdfixstarssel: {display name: True/False}
selmap = getattr(options, 'pdfixstarssel', None)
if isinstance(selmap, dict) and len(selmap) > 0:
name2code = _cat_name_to_code_map()
out = []
for nm, flag in selmap.items():
if not flag:
continue
nmu = nm if isinstance(nm, unicode) else unicode(nm)
code = name2code.get(nmu)
if code:
out.append(code)
if out:
return out
# 3) 아무 것도 못 찾음
return []
def _options_selected_stars(options):
"""
사용자가 켠 항성만 반환.
반환값: ['alTau', 'alLeo', ...] # fixstars.cat의 'code' 필드 (콤마 없이)
(paranwnd.py와 동일한 방식: options.fixstars 딕셔너리의 key만 사용)
"""
try:
fs = getattr(options, 'fixstars', {})
except Exception:
return []
# Py2/3 호환 키 이터레이터
try:
keys_iter = fs.iterkeys()
except AttributeError:
keys_iter = fs.keys()
out = []
for k in keys_iter:
# unicode 보정
try:
ku = k if isinstance(k, unicode) else unicode(k)
except NameError:
ku = k
out.append(ku.strip())
return out
# 2-2) 딕셔너리 {이름: bool}
if isinstance(cand, dict):
return [ (k if isinstance(k, unicode) else unicode(k)) for k, v in cand.items() if bool(v) ]
# 2-3) 리스트/튜플: (이름,) 또는 (이름, bool) 또는 그냥 이름 문자열
if isinstance(cand, (list, tuple)):
out = []
for it in cand:
if isinstance(it, (list, tuple)) and len(it) >= 1:
nm = it[0]
flag = True if len(it) == 1 else bool(it[1])
if flag:
out.append(nm if isinstance(nm, unicode) else unicode(nm))
elif isinstance(it, basestring):
out.append(it if isinstance(it, unicode) else unicode(it))
return out
# 2-4) 없거나 인식 실패 → 빈 목록
return []
# 3-2) 딕셔너리 {이름: bool}
if isinstance(cand, dict):
out = [ (k if isinstance(k, unicode) else unicode(k)) for k,v in cand.items() if v ]
return out
# 3-3) 리스트/튜플: (이름,) 또는 (이름, bool) 또는 그냥 이름 문자열들
if isinstance(cand, (list, tuple)):
out = []
for it in cand:
if isinstance(it, (list, tuple)) and len(it) >= 1:
nm = it[0]
flag = True if len(it) == 1 else bool(it[1])
if flag:
out.append(nm if isinstance(nm, unicode) else unicode(nm))
elif isinstance(it, basestring):
out.append(it if isinstance(it, unicode) else unicode(it))
return out
# 그 외는 지원 안 함 → 빈 목록
return []
def _fallback_names_from_cat(path=None):
"""
선택 리스트를 못 얻었을 때 최후의 안전장치.
Morinus의 ephe 경로에 있는 fixedstars.cat에서 '이름'만 긁어온다.
라인 포맷이 버전별로 다르므로 세미콜론/탭/다중 공백을 모두 시도.
"""
if path is None:
path = _fixstars_cat_path()
if not path:
return []
names = []
try:
with open(path, "r") as f:
for ln in f:
ln = ln.strip()
if not ln or ln.startswith("#"):
continue
if ";" in ln:
parts = [x.strip() for x in ln.split(";")]
nm = parts[0]
elif "\t" in ln:
nm = ln.split("\t", 1)[0].strip()
else:
# 연속 공백 분리
toks = ln.split()
nm = toks[0]
if nm:
names.append(unicode(nm) if not isinstance(nm, unicode) else nm)
except Exception:
pass
return names
def _star_ofdate_ra_dec(star_name, jd_ut):
"""
카탈로그(콤마/세미콜론/공백)에서 좌표를 읽고,
B1950이면 FK4→FK5(J2000) 변환 후 of-date 프리세션,
J2000/ICRS이면 바로 of-date 프리세션.
반환: of-date RA[deg], Dec[deg], "Name, mag"
"""
ra_cat, dec_cat, mag, frame = _cat_lookup_equ_generic(star_name)
if frame == 'B1950':
ra2000, dec2000 = _fk4_b1950_to_fk5_j2000(ra_cat, dec_cat)
else:
ra2000, dec2000 = ra_cat, dec_cat
# TT≈UT 가정(전통 PD 해상도엔 충분)
jd_tt = jd_ut + swe.swe_deltat(jd_ut)
ra_d, dec_d = _j2000_to_ofdate(ra2000, dec2000, jd_tt)
info = u"{0}, {1}".format(star_name, u"" if mag is None else u"{0}".format(mag)).strip()
return ra_d, dec_d, info
def _cat_lookup_equ_generic(star_name, path=None):
"""
fixstars/fixedstars.cat 의 좌표 파싱 (콤마/세미콜론/공백 포맷 모두 지원).
반환: (ra_deg, dec_deg, mag, frame) # frame: 'B1950' 또는 'J2000'
- 콤마 포맷 예:
Aldebaran ,alTau,ICRS,04,35,55.2387,16,30,33.485, ...
(name, alias, frame, RA_h,RA_m,RA_s, Dec_d,Dec_m,Dec_s, ...)
- 세미콜론 포맷 예:
Aldebaran; 04 35 55.24; +16 30 33.5; 0.86
- 공백/탭 포맷 예:
Aldebaran 04 35 55.24 +16 30 33.5 0.86
"""
if path is None:
path = _fixstars_cat_path()
if not path:
raise ValueError("fixed stars catalog not found")
def hms_to_deg(h, m, s): return (h + m/60.0 + s/3600.0) * 15.0
def dms_to_deg(sign, d, m, s):
val = d + m/60.0 + s/3600.0
return val if sign >= 0 else -val
want = star_name.strip().lower()
f = open(path, "r")
try:
for ln in f:
if not isinstance(ln, unicode):
try: ln = ln.decode('utf-8', 'ignore')
except: pass
s = ln.strip()
if (not s) or s.startswith("#"):
continue
name = None; alias = None; ra_tokens = None; dec_tokens = None
mag = None; frame = 'J2000' # 기본값(J2000/ICRS 취급)
if "," in s:
parts = [x.strip() for x in s.split(",")]
if len(parts) >= 9:
name = parts[0]
alias = parts[1] if len(parts) > 1 else None
fr = parts[2].upper() if len(parts) > 2 else ""
if "1950" in fr or "B1950" in fr or "FK4" in fr:
frame = 'B1950'
else:
frame = 'J2000' # 'ICRS','2000' 등은 전부 J2000 취급
ra_tokens = [parts[3], parts[4], parts[5]]
dec_tokens = [parts[6], parts[7], parts[8]]
# 뒤쪽 토큰들 중 -2~+7 사이 값(대략 Vmag) 하나를 잡아본다(선택)
for tok in parts[9:]:
tok = tok.strip()
try:
v = float(tok)
if -2.5 <= v <= 8.5:
mag = v
break
except:
pass
elif ";" in s:
parts = [x.strip() for x in s.split(";")]
if len(parts) >= 3:
name = parts[0]
ra_tokens = parts[1].replace(":", " ").split()
dec_tokens = parts[2].replace(":", " ").split()
frame = 'B1950'
if len(parts) >= 4:
try: mag = float(parts[3].split()[0])
except: pass
else:
toks = re.split(r"[ \t]+", s)
if len(toks) >= 7:
name = toks[0]
ra_tokens = [toks[1], toks[2], toks[3]]
dec_tokens = [toks[4], toks[5], toks[6]]
frame = 'B1950'
if len(toks) >= 8:
try: mag = float(toks[7])
except: pass
if not name or not ra_tokens or not dec_tokens:
continue
# 이름/별칭 매칭(대소문자/공백 무시)
nm = name.strip().lower()
al = alias.strip().lower() if alias else None
if want != nm and (al is None or want != al):
continue
# 숫자화
try:
rh, rm, rs = float(ra_tokens[0]), float(ra_tokens[1]), float(ra_tokens[2])
dd0 = dec_tokens[0]; sgn = -1 if dd0.startswith("-") else 1
dd = float(dd0.replace("+","").replace("-",""))
dm, ds = float(dec_tokens[1]), float(dec_tokens[2])
except:
continue
ra_deg = hms_to_deg(rh, rm, rs)
dec_deg = dms_to_deg(sgn, dd, dm, ds)
return ra_deg, dec_deg, mag, frame
finally:
try: f.close()
except: pass
raise ValueError("Star not found in catalog: {0}".format(star_name))
def _date_string_from_jd(jd_ut, chrt, options):
gregflag = _calendar_flag(chrt, options) # 1=Greg, 0=Jul
y, m, d, frac = swe.swe_revjul(jd_ut, gregflag) # UTC 기준
hh = int(frac * 24.0)
mm = int((frac * 24.0 - hh) * 60.0)
return u"{0:04d}-{1:02d}-{2:02d} {3:02d}:{4:02d} UTC".format(y, m, d, hh, mm)
def compute_fixedstar_angle_rows(horoscope, options, age_range=None, direction=None):
rows = []
# 0) 기본 상수
NAIBOD_DEG_PER_DAY = 360.0 / 365.2421897
# 1) 출생 RAMC/위도
ramc0 = _get_ramc0_deg(horoscope) # 기존 함수 그대로 사용
phi = _observer_lat(horoscope) # 위도 얻기 (이미 정의돼 있음)
jd0 = getattr(getattr(horoscope, 'time', None), 'jd', None)
if jd0 is None:
return []
if age_range is None:
age_range = (0.0, primdirs.PrimDirs.LIMIT)
if direction is None:
direction = primdirs.PrimDirs.BOTHDC
lo_age, hi_age = age_range
# 2) 사용자가 켠 항성 '코드'만 모으기
codes = _selected_star_codes(options)
if not codes:
return []
# 표시용 이름 가져오기 위해 DB 미리 로드
_db_cache = _load_fixstars_cat()
# 3) 각 항성(코드)에 대해 네 각도 계산 → 150년 확장
max_days = hi_age * 365.2421897
for code in codes:
# of-date RA/Dec 계산(세차+고유운동 반영)
try:
ra, dec, dispname = _ra_dec_star_ofdate_from_code(code, jd0)
except Exception:
continue
# 표시용 이름이 DB에 있으면 그걸로 덮어쓰기
try:
dispname = _db_cache.get(code, {}).get('name', dispname)
except Exception:
pass
if not isinstance(dispname, unicode):
try: dispname = unicode(dispname)
except: pass
# 옵션이 비어있으면 JSON에서 복구
if (not hasattr(options, 'fixstarAliasMap')) or (not isinstance(options.fixstarAliasMap, dict)) or (not options.fixstarAliasMap):
import common, os, json
alias_json = os.path.join(common.common.ephepath, 'fixstar_aliases.json')
if os.path.isfile(alias_json):
with open(alias_json, 'r') as _f:
_data = json.load(_f)
if isinstance(_data, dict):
options.fixstarAliasMap = dict(_data)
# [ADD] FixStarsDlg에서 사용자가 고른 표기(코드→표시명 맵)가 있으면 그것을 우선
try:
if hasattr(options, 'fixstarAliasMap') and code in options.fixstarAliasMap:
dispname = options.fixstarAliasMap[code]
except Exception:
pass
# 출생 RAMC 세트 (MC, IC, AOASC, DODESC)
ramc, raic, aoasc, dodesc = _ramc_pack(horoscope)
# ADlat로 ASC/DSC 가능성 판단
ad = _adlat(phi, dec)
# 선택된 프라이머리 디렉션 방식 (UTP는 항성-앵글에선 배제 → Placidian(Semiarc)로 취급)
pd_method = getattr(options, 'primarydir', primdirs.PrimDirs.PLACIDIANSEMIARC)
if pd_method == primdirs.PrimDirs.PLACIDIANUNDERTHEPOLE:
pd_method = primdirs.PrimDirs.PLACIDIANSEMIARC
# 각도쌍: (표시각, target, base)
# ==> 프라이머리의 옵션(Options.sigangles: [Asc, Dsc, MC, IC])을 그대로 적용
sigflags = getattr(options, 'sigangles', [True, True, True, True])
want = set()
try:
# PD와 동일 인덱스: [Asc, Dsc, MC, IC] = [0,1,2,3]
if sigflags[2]: want.add("MC")
if sigflags[3]: want.add("IC")
if sigflags[0]: want.add("ASC")
if sigflags[1]: want.add("DSC")
except Exception:
# 안전망: 모두 표시
want = {"MC", "IC", "ASC", "DSC"}
pairs = []
# MC/IC는 언제나 정의됨
if "MC" in want:
pairs.append(("MC", ra, ramc))
if "IC" in want:
pairs.append(("IC", ra, raic))
# ASC/DSC는 ADlat이 정의될 때만 가능
if ad is not None:
aostar = (ra - ad) % 360.0
dostar = (ra + ad) % 360.0
if "ASC" in want:
pairs.append(("ASC", aostar, aoasc))
if "DSC" in want:
pairs.append(("DSC", dostar, dodesc))
# 정방/역방 각각 1개씩만 생성 (k회전 확장 금지)
for sig, target, base in pairs:
arcD = _arc_direct(base, target) # 0..360
arcC = _arc_converse(base, target) # 0..360
# 이 두 줄을 아래 두 줄로 교체
yrsD = _arc_to_years_from_primary_key(horoscope, options, arcD, True)
yrsC = _arc_to_years_from_primary_key(horoscope, options, arcC, False)
# 선택한 방향/나이 범위에 맞춰 D/C 개별 추가
if direction in (primdirs.PrimDirs.BOTHDC, primdirs.PrimDirs.DIRECT):
if lo_age <= yrsD <= hi_age:
jd_evt = jd0 + yrsD * 365.2421897
rows.append({
'prom': dispname,
'dc' : 'D',
'sig' : sig,
'arc' : arcD,
'jd' : jd_evt,
})
if direction in (primdirs.PrimDirs.BOTHDC, primdirs.PrimDirs.CONVERSE):
if lo_age <= yrsC <= hi_age:
jd_evt = jd0 + yrsC * 365.2421897
rows.append({
'prom': dispname,
'dc' : 'C',
'sig' : sig,
'arc' : arcC,
'jd' : jd_evt,
})
# 4) 날짜 오름차순 정렬 후 문자열화
rows.sort(key=lambda r: r['jd'])
for r in rows:
r['date'] = _date_string_from_jd(r['jd'], horoscope, options)
return rows