The requested changes are too big and content was truncated. Show full diff
This diff has been collapsed as it changes many lines, (1420 lines changed) Show them Hide them | |||||
@@ -0,0 +1,1420 | |||||
|
1 | """ | |||
|
2 | The module ASTRO_COORDS.py gathers classes and functions for coordinates transformation. Additiona- | |||
|
3 | lly a class EquatorialCorrections and celestial bodies are defined. The first of these is to correct | |||
|
4 | any error in the location of the body and the second to know the location of certain celestial bo- | |||
|
5 | dies in the sky. | |||
|
6 | ||||
|
7 | MODULES CALLED: | |||
|
8 | OS, NUMPY, NUMERIC, SCIPY, TIME_CONVERSIONS | |||
|
9 | ||||
|
10 | MODIFICATION HISTORY: | |||
|
11 | Created by Ing. Freddy Galindo (frederickgalindo@gmail.com). ROJ Sep 20, 2009. | |||
|
12 | """ | |||
|
13 | ||||
|
14 | import numpy | |||
|
15 | #import Numeric | |||
|
16 | import scipy.interpolate | |||
|
17 | import os | |||
|
18 | import sys | |||
|
19 | from utils import TimeTools | |||
|
20 | from utils import Misc_Routines | |||
|
21 | ||||
|
22 | class EquatorialCorrections(): | |||
|
23 | def __init__(self): | |||
|
24 | """ | |||
|
25 | EquatorialCorrections class creates an object to call methods to correct the loca- | |||
|
26 | tion of the celestial bodies. | |||
|
27 | ||||
|
28 | Modification History | |||
|
29 | -------------------- | |||
|
30 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 27 September 2009. | |||
|
31 | """ | |||
|
32 | ||||
|
33 | pass | |||
|
34 | ||||
|
35 | def co_nutate(self,jd,ra,dec): | |||
|
36 | """ | |||
|
37 | co_nutate calculates changes in RA and Dec due to nutation of the Earth's rotation | |||
|
38 | Additionally it returns the obliquity of the ecliptic (eps), nutation in the longi- | |||
|
39 | tude of the ecliptic (d_psi) and nutation in the pbliquity of the ecliptic (d_eps). | |||
|
40 | ||||
|
41 | Parameters | |||
|
42 | ---------- | |||
|
43 | jd = Julian Date (Scalar or array). | |||
|
44 | RA = A scalar o array giving the Right Ascention of interest. | |||
|
45 | Dec = A scalar o array giving the Right Ascention of interest. | |||
|
46 | ||||
|
47 | Return | |||
|
48 | ------ | |||
|
49 | d_ra = Correction to ra due to nutation. | |||
|
50 | d_dec = Correction to dec due to nutation. | |||
|
51 | ||||
|
52 | Examples | |||
|
53 | -------- | |||
|
54 | >> Julian = 2462088.7 | |||
|
55 | >> Ra = 41.547213 | |||
|
56 | >> Dec = 49.348483 | |||
|
57 | >> [d_ra,d_dec,eps,d_psi,d_eps] = co_nutate(julian,Ra,Dec) | |||
|
58 | >> print d_ra, d_dec, eps, d_psi, d_eps | |||
|
59 | [ 15.84276651] [ 6.21641029] [ 0.4090404] [ 14.85990198] [ 2.70408658] | |||
|
60 | ||||
|
61 | Modification history | |||
|
62 | -------------------- | |||
|
63 | Written by Chris O'Dell, 2002. | |||
|
64 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
65 | """ | |||
|
66 | ||||
|
67 | jd = numpy.atleast_1d(jd) | |||
|
68 | ra = numpy.atleast_1d(ra) | |||
|
69 | dec = numpy.atleast_1d(dec) | |||
|
70 | ||||
|
71 | # Useful transformation constants | |||
|
72 | d2as = numpy.pi/(180.*3600.) | |||
|
73 | ||||
|
74 | # Julian centuries from J2000 of jd | |||
|
75 | T = (jd - 2451545.0)/36525.0 | |||
|
76 | ||||
|
77 | # Must calculate obliquity of ecliptic | |||
|
78 | [d_psi, d_eps] = self.nutate(jd) | |||
|
79 | d_psi = numpy.atleast_1d(d_psi) | |||
|
80 | d_eps = numpy.atleast_1d(d_eps) | |||
|
81 | ||||
|
82 | eps0 = (23.4392911*3600.) - (46.8150*T) - (0.00059*T**2) + (0.001813*T**3) | |||
|
83 | # True obliquity of the ecliptic in radians | |||
|
84 | eps = (eps0 + d_eps)/3600.*Misc_Routines.CoFactors.d2r | |||
|
85 | ||||
|
86 | # Useful numbers | |||
|
87 | ce = numpy.cos(eps) | |||
|
88 | se = numpy.sin(eps) | |||
|
89 | ||||
|
90 | # Convert Ra-Dec to equatorial rectangular coordinates | |||
|
91 | x = numpy.cos(ra*Misc_Routines.CoFactors.d2r)*numpy.cos(dec*Misc_Routines.CoFactors.d2r) | |||
|
92 | y = numpy.sin(ra*Misc_Routines.CoFactors.d2r)*numpy.cos(dec*Misc_Routines.CoFactors.d2r) | |||
|
93 | z = numpy.sin(dec*Misc_Routines.CoFactors.d2r) | |||
|
94 | ||||
|
95 | # Apply corrections to each rectangular coordinate | |||
|
96 | x2 = x - (y*ce + z*se)*d_psi*Misc_Routines.CoFactors.s2r | |||
|
97 | y2 = y + (x*ce*d_psi - z*d_eps)*Misc_Routines.CoFactors.s2r | |||
|
98 | z2 = z + (x*se*d_psi + y*d_eps)*Misc_Routines.CoFactors.s2r | |||
|
99 | ||||
|
100 | # Convert bask to equatorial spherical coordinates | |||
|
101 | r = numpy.sqrt(x2**2. + y2**2. + z2**2.) | |||
|
102 | xyproj =numpy.sqrt(x2**2. + y2**2.) | |||
|
103 | ||||
|
104 | ra2 = x2*0.0 | |||
|
105 | dec2 = x2*0.0 | |||
|
106 | ||||
|
107 | xyproj = numpy.atleast_1d(xyproj) | |||
|
108 | z = numpy.atleast_1d(z) | |||
|
109 | r = numpy.atleast_1d(r) | |||
|
110 | x2 = numpy.atleast_1d(x2) | |||
|
111 | y2 = numpy.atleast_1d(y2) | |||
|
112 | z2 = numpy.atleast_1d(z2) | |||
|
113 | ra2 = numpy.atleast_1d(ra2) | |||
|
114 | dec2 = numpy.atleast_1d(dec2) | |||
|
115 | ||||
|
116 | w1 = numpy.where((xyproj==0) & (z!=0)) | |||
|
117 | w2 = numpy.where(xyproj!=0) | |||
|
118 | ||||
|
119 | # Calculate Ra and Dec in radians (later convert to degrees) | |||
|
120 | if w1[0].size>0: | |||
|
121 | # Places where xyproj=0 (point at NCP or SCP) | |||
|
122 | dec2[w1] = numpy.arcsin(z2[w1]/r[w1]) | |||
|
123 | ra2[w1] = 0 | |||
|
124 | ||||
|
125 | if w2[0].size>0: | |||
|
126 | # Places other than NCP or SCP | |||
|
127 | ra2[w2] = numpy.arctan2(y2[w2],x2[w2]) | |||
|
128 | dec2[w2] = numpy.arcsin(z2[w2]/r[w2]) | |||
|
129 | ||||
|
130 | # Converting to degree | |||
|
131 | ra2 = ra2/Misc_Routines.CoFactors.d2r | |||
|
132 | dec2 = dec2/Misc_Routines.CoFactors.d2r | |||
|
133 | ||||
|
134 | w = numpy.where(ra2<0.) | |||
|
135 | if w[0].size>0: | |||
|
136 | ra2[w] = ra2[w] + 360. | |||
|
137 | ||||
|
138 | # Return changes in Ra and Dec in arcseconds | |||
|
139 | d_ra = (ra2 -ra)*3600. | |||
|
140 | d_dec = (dec2 - dec)*3600. | |||
|
141 | ||||
|
142 | return d_ra, d_dec, eps, d_psi, d_eps | |||
|
143 | ||||
|
144 | def nutate(self,jd): | |||
|
145 | """ | |||
|
146 | nutate returns the nutation in longitude and obliquity for a given Julian date. | |||
|
147 | ||||
|
148 | Parameters | |||
|
149 | ---------- | |||
|
150 | jd = Julian ephemeris date, scalar or vector. | |||
|
151 | ||||
|
152 | Return | |||
|
153 | ------ | |||
|
154 | nut_long = The nutation in longitude. | |||
|
155 | nut_obliq = The nutation in latitude. | |||
|
156 | ||||
|
157 | Example | |||
|
158 | ------- | |||
|
159 | >> julian = 2446895.5 | |||
|
160 | >> [nut_long,nut_obliq] = nutate(julian) | |||
|
161 | >> print nut_long, nut_obliq | |||
|
162 | -3.78793107711 9.44252069864 | |||
|
163 | ||||
|
164 | >> julians = 2415020.5 + numpy.arange(50) | |||
|
165 | >> [nut_long,nut_obliq] = nutate(julians) | |||
|
166 | ||||
|
167 | Modification History | |||
|
168 | -------------------- | |||
|
169 | Written by W.Landsman (Goddard/HSTX), June 1996. | |||
|
170 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
171 | """ | |||
|
172 | ||||
|
173 | jd = numpy.atleast_1d(jd) | |||
|
174 | ||||
|
175 | # Form time in Julian centuries from 1900 | |||
|
176 | t = (jd - 2451545.0)/36525.0 | |||
|
177 | ||||
|
178 | # Mean elongation of the moon | |||
|
179 | coeff1 = numpy.array([1/189474.0,-0.0019142,445267.111480,297.85036]) | |||
|
180 | d = numpy.poly1d(coeff1) | |||
|
181 | d = d(t)*Misc_Routines.CoFactors.d2r | |||
|
182 | d = self.cirrange(d,rad=1) | |||
|
183 | ||||
|
184 | # Sun's mean elongation | |||
|
185 | coeff2 = numpy.array([-1./3e5,-0.0001603,35999.050340,357.52772]) | |||
|
186 | m = numpy.poly1d(coeff2) | |||
|
187 | m = m(t)*Misc_Routines.CoFactors.d2r | |||
|
188 | m = self.cirrange(m,rad=1) | |||
|
189 | ||||
|
190 | # Moon's mean elongation | |||
|
191 | coeff3 = numpy.array([1.0/5.625e4,0.0086972,477198.867398,134.96298]) | |||
|
192 | mprime = numpy.poly1d(coeff3) | |||
|
193 | mprime = mprime(t)*Misc_Routines.CoFactors.d2r | |||
|
194 | mprime = self.cirrange(mprime,rad=1) | |||
|
195 | ||||
|
196 | # Moon's argument of latitude | |||
|
197 | coeff4 = numpy.array([-1.0/3.27270e5,-0.0036825,483202.017538,93.27191]) | |||
|
198 | f = numpy.poly1d(coeff4) | |||
|
199 | f = f(t)*Misc_Routines.CoFactors.d2r | |||
|
200 | f = self.cirrange(f,rad=1) | |||
|
201 | ||||
|
202 | # Longitude fo the ascending node of the Moon's mean orbit on the ecliptic, measu- | |||
|
203 | # red from the mean equinox of the date. | |||
|
204 | coeff5 = numpy.array([1.0/4.5e5,0.0020708,-1934.136261,125.04452]) | |||
|
205 | omega = numpy.poly1d(coeff5) | |||
|
206 | omega = omega(t)*Misc_Routines.CoFactors.d2r | |||
|
207 | omega = self.cirrange(omega,rad=1) | |||
|
208 | ||||
|
209 | d_lng = numpy.array([0,-2,0,0,0,0,-2,0,0,-2,-2,-2,0,2,0,2,0,0,-2,0,2,0,0,-2,0,-2,0,0,\ | |||
|
210 | 2,-2,0,-2,0,0,2,2,0,-2,0,2,2,-2,-2,2,2,0,-2,-2,0,-2,-2,0,-1,-2,1,0,0,-1,0,\ | |||
|
211 | 0,2,0,2]) | |||
|
212 | ||||
|
213 | m_lng = numpy.array([0,0,0,0,1,0,1,0,0,-1]) | |||
|
214 | m_lng = numpy.append(m_lng,numpy.zeros(17)) | |||
|
215 | m_lng = numpy.append(m_lng,numpy.array([2,0,2,1,0,-1,0,0,0,1,1,-1,0,0,0,0,0,0,-1,-1,0,0,\ | |||
|
216 | 0,1,0,0,1,0,0,0,-1,1,-1,-1,0,-1])) | |||
|
217 | ||||
|
218 | mp_lng = numpy.array([0,0,0,0,0,1,0,0,1,0,1,0,-1,0,1,-1,-1,1,2,-2,0,2,2,1,0,0, -1, 0,\ | |||
|
219 | -1,0,0,1,0,2,-1,1,0,1,0,0,1,2,1,-2,0,1,0,0,2,2,0,1,1,0,0,1,-2,1,1,1,-1,3,0]) | |||
|
220 | ||||
|
221 | f_lng = numpy.array([0,2,2,0,0,0,2,2,2,2,0,2,2,0,0,2,0,2,0,2,2,2,0,2,2,2,2,0,0,2,0,0,\ | |||
|
222 | 0,-2,2,2,2,0,2,2,0,2,2,0,0,0,2,0,2,0,2,-2,0,0,0,2,2,0,0,2,2,2,2]) | |||
|
223 | ||||
|
224 | om_lng = numpy.array([1,2,2,2,0,0,2,1,2,2,0,1,2,0,1,2,1,1,0,1,2,2,0,2,0,0,1,0,1,2,1, \ | |||
|
225 | 1,1,0,1,2,2,0,2,1,0,2,1,1,1,0,1,1,1,1,1,0,0,0,0,0,2,0,0,2,2,2,2]) | |||
|
226 | ||||
|
227 | sin_lng = numpy.array([-171996,-13187,-2274,2062,1426,712,-517,-386,-301, 217, -158, \ | |||
|
228 | 129,123,63,63,-59,-58,-51,48,46,-38,-31,29,29,26,-22,21,17,16,-16,-15,-13,\ | |||
|
229 | -12,11,-10,-8,7,-7,-7,-7,6,6,6,-6,-6,5,-5,-5,-5,4,4,4,-4,-4,-4,3,-3,-3,-3,\ | |||
|
230 | -3,-3,-3,-3]) | |||
|
231 | ||||
|
232 | sdelt = numpy.array([-174.2,-1.6,-0.2,0.2,-3.4,0.1,1.2,-0.4,0,-0.5,0, 0.1, 0, 0, 0.1,\ | |||
|
233 | 0,-0.1]) | |||
|
234 | sdelt = numpy.append(sdelt,numpy.zeros(10)) | |||
|
235 | sdelt = numpy.append(sdelt,numpy.array([-0.1, 0, 0.1])) | |||
|
236 | sdelt = numpy.append(sdelt,numpy.zeros(33)) | |||
|
237 | ||||
|
238 | cos_lng = numpy.array([92025,5736,977,-895,54,-7,224,200,129,-95,0,-70,-53,0,-33,26, \ | |||
|
239 | 32,27,0,-24,16,13,0,-12,0,0,-10,0,-8,7,9,7,6,0,5,3,-3,0,3,3,0,-3,-3,3,3,0,\ | |||
|
240 | 3,3,3]) | |||
|
241 | cos_lng = numpy.append(cos_lng,numpy.zeros(14)) | |||
|
242 | ||||
|
243 | cdelt = numpy.array([8.9,-3.1,-0.5,0.5,-0.1,0.0,-0.6,0.0,-0.1,0.3]) | |||
|
244 | cdelt = numpy.append(cdelt,numpy.zeros(53)) | |||
|
245 | ||||
|
246 | # Sum the periodic terms. | |||
|
247 | n = numpy.size(jd) | |||
|
248 | nut_long = numpy.zeros(n) | |||
|
249 | nut_obliq = numpy.zeros(n) | |||
|
250 | ||||
|
251 | d_lng = d_lng.reshape(numpy.size(d_lng),1) | |||
|
252 | d = d.reshape(numpy.size(d),1) | |||
|
253 | matrix_d_lng = numpy.dot(d_lng,d.transpose()) | |||
|
254 | ||||
|
255 | m_lng = m_lng.reshape(numpy.size(m_lng),1) | |||
|
256 | m = m.reshape(numpy.size(m),1) | |||
|
257 | matrix_m_lng = numpy.dot(m_lng,m.transpose()) | |||
|
258 | ||||
|
259 | mp_lng = mp_lng.reshape(numpy.size(mp_lng),1) | |||
|
260 | mprime = mprime.reshape(numpy.size(mprime),1) | |||
|
261 | matrix_mp_lng = numpy.dot(mp_lng,mprime.transpose()) | |||
|
262 | ||||
|
263 | f_lng = f_lng.reshape(numpy.size(f_lng),1) | |||
|
264 | f = f.reshape(numpy.size(f),1) | |||
|
265 | matrix_f_lng = numpy.dot(f_lng,f.transpose()) | |||
|
266 | ||||
|
267 | om_lng = om_lng.reshape(numpy.size(om_lng),1) | |||
|
268 | omega = omega.reshape(numpy.size(omega),1) | |||
|
269 | matrix_om_lng = numpy.dot(om_lng,omega.transpose()) | |||
|
270 | ||||
|
271 | arg = matrix_d_lng + matrix_m_lng + matrix_mp_lng + matrix_f_lng + matrix_om_lng | |||
|
272 | ||||
|
273 | sarg = numpy.sin(arg) | |||
|
274 | carg = numpy.cos(arg) | |||
|
275 | ||||
|
276 | for ii in numpy.arange(n): | |||
|
277 | nut_long[ii] = 0.0001*numpy.sum((sdelt*t[ii] + sin_lng)*sarg[:,ii]) | |||
|
278 | nut_obliq[ii] = 0.0001*numpy.sum((cdelt*t[ii] + cos_lng)*carg[:,ii]) | |||
|
279 | ||||
|
280 | if numpy.size(jd)==1: | |||
|
281 | nut_long = nut_long[0] | |||
|
282 | nut_obliq = nut_obliq[0] | |||
|
283 | ||||
|
284 | return nut_long, nut_obliq | |||
|
285 | ||||
|
286 | def co_aberration(self,jd,ra,dec): | |||
|
287 | """ | |||
|
288 | co_aberration calculates changes to Ra and Dec due to "the effect of aberration". | |||
|
289 | ||||
|
290 | Parameters | |||
|
291 | ---------- | |||
|
292 | jd = Julian Date (Scalar or vector). | |||
|
293 | ra = A scalar o vector giving the Right Ascention of interest. | |||
|
294 | dec = A scalar o vector giving the Declination of interest. | |||
|
295 | ||||
|
296 | Return | |||
|
297 | ------ | |||
|
298 | d_ra = The correction to right ascension due to aberration (must be added to ra to | |||
|
299 | get the correct value). | |||
|
300 | d_dec = The correction to declination due to aberration (must be added to the dec | |||
|
301 | to get the correct value). | |||
|
302 | eps = True obliquity of the ecliptic (in radians). | |||
|
303 | ||||
|
304 | Examples | |||
|
305 | -------- | |||
|
306 | >> Julian = 2462088.7 | |||
|
307 | >> Ra = 41.547213 | |||
|
308 | >> Dec = 49.348483 | |||
|
309 | >> [d_ra,d_dec,eps] = co_aberration(julian,Ra,Dec) | |||
|
310 | >> print d_ra, d_dec, eps | |||
|
311 | [ 30.04441796] [ 6.69837858] [ 0.40904059] | |||
|
312 | ||||
|
313 | Modification history | |||
|
314 | -------------------- | |||
|
315 | Written by Chris O'Dell , Univ. of Wisconsin, June 2002. | |||
|
316 | Converted to Python by Freddy R. Galindo, ROJ, 27 September 2009. | |||
|
317 | """ | |||
|
318 | ||||
|
319 | # Julian centuries from J2000 of jd. | |||
|
320 | T = (jd - 2451545.0)/36525.0 | |||
|
321 | ||||
|
322 | # Getting obliquity of ecliptic | |||
|
323 | njd = numpy.size(jd) | |||
|
324 | jd = numpy.atleast_1d(jd) | |||
|
325 | ra = numpy.atleast_1d(ra) | |||
|
326 | dec = numpy.atleast_1d(dec) | |||
|
327 | ||||
|
328 | d_psi = numpy.zeros(njd) | |||
|
329 | d_epsilon = d_psi | |||
|
330 | for ii in numpy.arange(njd): | |||
|
331 | [dp,de] = self.nutate(jd[ii]) | |||
|
332 | d_psi[ii] = dp | |||
|
333 | d_epsilon[ii] = de | |||
|
334 | ||||
|
335 | coeff = 23 + 26/60. + 21.488/3600. | |||
|
336 | eps0 = coeff*3600. - 46.8150*T - 0.00059*T**2. + 0.001813*T**3. | |||
|
337 | # True obliquity of the ecliptic in radians | |||
|
338 | eps = (eps0 + d_epsilon)/3600*Misc_Routines.CoFactors.d2r | |||
|
339 | ||||
|
340 | celestialbodies = CelestialBodies() | |||
|
341 | [sunra,sundec,sunlon,sunobliq] = celestialbodies.sunpos(jd) | |||
|
342 | ||||
|
343 | # Earth's orbital eccentricity | |||
|
344 | e = 0.016708634 - 0.000042037*T - 0.0000001267*T**2. | |||
|
345 | ||||
|
346 | # longitude of perihelion, in degrees | |||
|
347 | pi = 102.93735 + 1.71946*T + 0.00046*T**2 | |||
|
348 | ||||
|
349 | # Constant of aberration, in arcseconds | |||
|
350 | k = 20.49552 | |||
|
351 | ||||
|
352 | cd = numpy.cos(dec*Misc_Routines.CoFactors.d2r) ; sd = numpy.sin(dec*Misc_Routines.CoFactors.d2r) | |||
|
353 | ce = numpy.cos(eps) ; te = numpy.tan(eps) | |||
|
354 | cp = numpy.cos(pi*Misc_Routines.CoFactors.d2r) ; sp = numpy.sin(pi*Misc_Routines.CoFactors.d2r) | |||
|
355 | cs = numpy.cos(sunlon*Misc_Routines.CoFactors.d2r) ; ss = numpy.sin(sunlon*Misc_Routines.CoFactors.d2r) | |||
|
356 | ca = numpy.cos(ra*Misc_Routines.CoFactors.d2r) ; sa = numpy.sin(ra*Misc_Routines.CoFactors.d2r) | |||
|
357 | ||||
|
358 | term1 = (ca*cs*ce + sa*ss)/cd | |||
|
359 | term2 = (ca*cp*ce + sa*sp)/cd | |||
|
360 | term3 = (cs*ce*(te*cd - sa*sd) + ca*sd*ss) | |||
|
361 | term4 = (cp*ce*(te*cd - sa*sd) + ca*sd*sp) | |||
|
362 | ||||
|
363 | d_ra = -k*term1 + e*k*term2 | |||
|
364 | d_dec = -k*term3 + e*k*term4 | |||
|
365 | ||||
|
366 | return d_ra, d_dec, eps | |||
|
367 | ||||
|
368 | def precess(self,ra,dec,equinox1=None,equinox2=None,FK4=0,rad=0): | |||
|
369 | """ | |||
|
370 | precess coordinates from EQUINOX1 to EQUINOX2 | |||
|
371 | ||||
|
372 | Parameters | |||
|
373 | ----------- | |||
|
374 | ra = A scalar o vector giving the Right Ascention of interest. | |||
|
375 | dec = A scalar o vector giving the Declination of interest. | |||
|
376 | equinox1 = Original equinox of coordinates, numeric scalar. If omitted, the __Pre- | |||
|
377 | cess will query for equinox1 and equinox2. | |||
|
378 | equinox2 = Original equinox of coordinates. | |||
|
379 | FK4 = If this keyword is set and non-zero, the FK4 (B1950) system will be used | |||
|
380 | otherwise FK5 (J2000) will be used instead. | |||
|
381 | rad = If this keyword is set and non-zero, then the input and output RAD and DEC | |||
|
382 | vectors are in radian rather than degree. | |||
|
383 | ||||
|
384 | Return | |||
|
385 | ------ | |||
|
386 | ra = Right ascension after precession (scalar or vector) in degrees, unless the rad | |||
|
387 | keyword is set. | |||
|
388 | dec = Declination after precession (scalar or vector) in degrees, unless the rad | |||
|
389 | keyword is set. | |||
|
390 | ||||
|
391 | Examples | |||
|
392 | -------- | |||
|
393 | >> Ra = 329.88772 | |||
|
394 | >> Dec = -56.992515 | |||
|
395 | >> [p_ra,p_dec] = precess(Ra,Dec,1950,1975,FK4=1) | |||
|
396 | >> print p_ra, p_dec | |||
|
397 | [ 330.31442971] [-56.87186154] | |||
|
398 | ||||
|
399 | Modification history | |||
|
400 | -------------------- | |||
|
401 | Written by Wayne Landsman, STI Corporation, August 1986. | |||
|
402 | Converted to Python by Freddy R. Galindo, ROJ, 27 September 2009. | |||
|
403 | """ | |||
|
404 | ||||
|
405 | npts = numpy.size(ra) | |||
|
406 | ra = numpy.atleast_1d(ra) | |||
|
407 | dec = numpy.atleast_1d(dec) | |||
|
408 | ||||
|
409 | if rad==0: | |||
|
410 | ra_rad = ra*Misc_Routines.CoFactors.d2r | |||
|
411 | dec_rad = dec*Misc_Routines.CoFactors.d2r | |||
|
412 | else: | |||
|
413 | ra_rad = ra | |||
|
414 | dec_rad = dec | |||
|
415 | ||||
|
416 | x = numpy.zeros((npts,3)) | |||
|
417 | x[:,0] = numpy.cos(dec_rad)*numpy.cos(ra_rad) | |||
|
418 | x[:,1] = numpy.cos(dec_rad)*numpy.sin(ra_rad) | |||
|
419 | x[:,2] = numpy.sin(dec_rad) | |||
|
420 | ||||
|
421 | # Use premat function to get precession matrix from equinox1 to equinox2 | |||
|
422 | r = self.premat(equinox1,equinox2,FK4) | |||
|
423 | ||||
|
424 | x2 = numpy.dot(r,x.transpose()) | |||
|
425 | ||||
|
426 | ra_rad = numpy.arctan2(x2[1,:],x2[0,:]) | |||
|
427 | dec_rad = numpy.arcsin(x2[2,:]) | |||
|
428 | ||||
|
429 | if rad==0: | |||
|
430 | ra = ra_rad/Misc_Routines.CoFactors.d2r | |||
|
431 | ra = ra + (ra<0)*360. | |||
|
432 | dec = dec_rad/Misc_Routines.CoFactors.d2r | |||
|
433 | else: | |||
|
434 | ra = ra_rad | |||
|
435 | ra = ra + (ra<0)*numpy.pi*2. | |||
|
436 | dec = dec_rad | |||
|
437 | ||||
|
438 | return ra, dec | |||
|
439 | ||||
|
440 | def premat(self,equinox1,equinox2,FK4=0): | |||
|
441 | """ | |||
|
442 | premat returns the precession matrix needed to go from EQUINOX1 to EQUINOX2. | |||
|
443 | ||||
|
444 | Parameters | |||
|
445 | ---------- | |||
|
446 | equinox1 = Original equinox of coordinates, numeric scalar. | |||
|
447 | equinox2 = Equinox of precessed coordinates. | |||
|
448 | FK4 = If this keyword is set and non-zero, the FK4 (B1950) system precession angles | |||
|
449 | are used to compute the precession matrix. The default is to use FK5 (J2000) pre- | |||
|
450 | cession angles. | |||
|
451 | ||||
|
452 | Return | |||
|
453 | ------ | |||
|
454 | r = Precession matrix, used to precess equatorial rectangular coordinates. | |||
|
455 | ||||
|
456 | Examples | |||
|
457 | -------- | |||
|
458 | >> matrix = premat(1950.0,1975.0,FK4=1) | |||
|
459 | >> print matrix | |||
|
460 | [[ 9.99981438e-01 -5.58774959e-03 -2.42908517e-03] | |||
|
461 | [ 5.58774959e-03 9.99984388e-01 -6.78691471e-06] | |||
|
462 | [ 2.42908517e-03 -6.78633095e-06 9.99997050e-01]] | |||
|
463 | ||||
|
464 | Modification history | |||
|
465 | -------------------- | |||
|
466 | Written by Wayne Landsman, HSTX Corporation, June 1994. | |||
|
467 | Converted to Python by Freddy R. Galindo, ROJ, 27 September 2009. | |||
|
468 | """ | |||
|
469 | ||||
|
470 | t = 0.001*(equinox2 - equinox1) | |||
|
471 | ||||
|
472 | if FK4==0: | |||
|
473 | st=0.001*(equinox1 - 2000.) | |||
|
474 | # Computing 3 rotation angles. | |||
|
475 | A=Misc_Routines.CoFactors.s2r*t*(23062.181+st*(139.656+0.0139*st)+t*(30.188-0.344*st+17.998*t)) | |||
|
476 | B=Misc_Routines.CoFactors.s2r*t*t*(79.280+0.410*st+0.205*t)+A | |||
|
477 | C=Misc_Routines.CoFactors.s2r*t*(20043.109-st*(85.33+0.217*st)+ t*(-42.665-0.217*st-41.833*t)) | |||
|
478 | else: | |||
|
479 | st=0.001*(equinox1 - 1900) | |||
|
480 | # Computing 3 rotation angles | |||
|
481 | A=Misc_Routines.CoFactors.s2r*t*(23042.53+st*(139.75+0.06*st)+t*(30.23-0.27*st+18.0*t)) | |||
|
482 | B=Misc_Routines.CoFactors.s2r*t*t*(79.27+0.66*st+0.32*t)+A | |||
|
483 | C=Misc_Routines.CoFactors.s2r*t*(20046.85-st*(85.33+0.37*st)+t*(-42.67-0.37*st-41.8*t)) | |||
|
484 | ||||
|
485 | sina = numpy.sin(A); sinb = numpy.sin(B); sinc = numpy.sin(C) | |||
|
486 | cosa = numpy.cos(A); cosb = numpy.cos(B); cosc = numpy.cos(C) | |||
|
487 | ||||
|
488 | r = numpy.zeros((3,3)) | |||
|
489 | r[:,0] = numpy.array([cosa*cosb*cosc-sina*sinb,sina*cosb+cosa*sinb*cosc,cosa*sinc]) | |||
|
490 | r[:,1] = numpy.array([-cosa*sinb-sina*cosb*cosc,cosa*cosb-sina*sinb*cosc,-sina*sinc]) | |||
|
491 | r[:,2] = numpy.array([-cosb*sinc,-sinb*sinc,cosc]) | |||
|
492 | ||||
|
493 | return r | |||
|
494 | ||||
|
495 | def cirrange(self,angle,rad=0): | |||
|
496 | """ | |||
|
497 | cirrange forces an angle into the range 0<= angle < 360. | |||
|
498 | ||||
|
499 | Parameters | |||
|
500 | ---------- | |||
|
501 | angle = The angle to modify, in degrees. Can be scalar or vector. | |||
|
502 | rad = Set to 1 if the angle is specified in radians rather than degrees. It is for- | |||
|
503 | ced into the range 0 <= angle < 2 PI | |||
|
504 | ||||
|
505 | Return | |||
|
506 | ------ | |||
|
507 | angle = The angle after the modification. | |||
|
508 | ||||
|
509 | Example | |||
|
510 | ------- | |||
|
511 | >> angle = cirrange(numpy.array([420,400,361])) | |||
|
512 | >> print angle | |||
|
513 | >> [60, 40, 1] | |||
|
514 | ||||
|
515 | Modification History | |||
|
516 | -------------------- | |||
|
517 | Written by Michael R. Greason, Hughes STX, 10 February 1994. | |||
|
518 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
519 | """ | |||
|
520 | ||||
|
521 | angle = numpy.atleast_1d(angle) | |||
|
522 | ||||
|
523 | if rad==1: | |||
|
524 | cnst = numpy.pi*2. | |||
|
525 | elif rad==0: | |||
|
526 | cnst = 360. | |||
|
527 | ||||
|
528 | # Deal with the lower limit. | |||
|
529 | angle = angle % cnst | |||
|
530 | ||||
|
531 | # Deal with negative values, if way | |||
|
532 | neg = numpy.where(angle<0.0) | |||
|
533 | if neg[0].size>0: angle[neg] = angle[neg] + cnst | |||
|
534 | ||||
|
535 | return angle | |||
|
536 | ||||
|
537 | ||||
|
538 | class CelestialBodies(EquatorialCorrections): | |||
|
539 | def __init__(self): | |||
|
540 | """ | |||
|
541 | CelestialBodies class creates a object to call methods of celestial bodies location. | |||
|
542 | ||||
|
543 | Modification History | |||
|
544 | -------------------- | |||
|
545 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 27 September 2009. | |||
|
546 | """ | |||
|
547 | ||||
|
548 | EquatorialCorrections.__init__(self) | |||
|
549 | ||||
|
550 | def sunpos(self,jd,rad=0): | |||
|
551 | """ | |||
|
552 | sunpos method computes the RA and Dec of the Sun at a given date. | |||
|
553 | ||||
|
554 | Parameters | |||
|
555 | ---------- | |||
|
556 | jd = The julian date of the day (and time), scalar or vector. | |||
|
557 | rad = If this keyword is set and non-zero, then the input and output RAD and DEC | |||
|
558 | vectors are in radian rather than degree. | |||
|
559 | ||||
|
560 | Return | |||
|
561 | ------ | |||
|
562 | ra = The right ascension of the sun at that date in degrees. | |||
|
563 | dec = The declination of the sun at that date in degrees. | |||
|
564 | elong = Ecliptic longitude of the sun at that date in degrees. | |||
|
565 | obliquity = The declination of the sun at that date in degrees. | |||
|
566 | ||||
|
567 | Examples | |||
|
568 | -------- | |||
|
569 | >> jd = 2466880 | |||
|
570 | >> [ra,dec,elong,obliquity] = sunpos(jd) | |||
|
571 | >> print ra, dec, elong, obliquity | |||
|
572 | [ 275.53499556] [-23.33840558] [ 275.08917968] [ 23.43596165] | |||
|
573 | ||||
|
574 | >> [ra,dec,elong,obliquity] = sunpos(jd,rad=1) | |||
|
575 | >> print ra, dec, elong, obliquity | |||
|
576 | [ 4.80899288] [-0.40733202] [ 4.80121192] [ 0.40903469] | |||
|
577 | ||||
|
578 | >> jd = 2450449.5 + numpy.arange(365) | |||
|
579 | >> [ra,dec,elong,obliquity] = sunpos(jd) | |||
|
580 | ||||
|
581 | Modification history | |||
|
582 | -------------------- | |||
|
583 | Written by Micheal R. Greason, STX Corporation, 28 October 1988. | |||
|
584 | Converted to Python by Freddy R. Galindo, ROJ, 27 September 2009. | |||
|
585 | """ | |||
|
586 | ||||
|
587 | jd = numpy.atleast_1d(jd) | |||
|
588 | ||||
|
589 | # Form time in Julian centuries from 1900. | |||
|
590 | t = (jd -2415020.0)/36525.0 | |||
|
591 | ||||
|
592 | # Form sun's mean longitude | |||
|
593 | l = (279.696678+((36000.768925*t) % 360.0))*3600.0 | |||
|
594 | ||||
|
595 | # Allow for ellipticity of the orbit (equation of centre) using the Earth's mean | |||
|
596 | # anomoly ME | |||
|
597 | me = 358.475844 + ((35999.049750*t) % 360.0) | |||
|
598 | ellcor = (6910.1 - 17.2*t)*numpy.sin(me*Misc_Routines.CoFactors.d2r) + 72.3*numpy.sin(2.0*me*Misc_Routines.CoFactors.d2r) | |||
|
599 | l = l + ellcor | |||
|
600 | ||||
|
601 | # Allow for the Venus perturbations using the mean anomaly of Venus MV | |||
|
602 | mv = 212.603219 + ((58517.803875*t) % 360.0) | |||
|
603 | vencorr = 4.8*numpy.cos((299.1017 + mv - me)*Misc_Routines.CoFactors.d2r) + \ | |||
|
604 | 5.5*numpy.cos((148.3133 + 2.0*mv - 2.0*me )*Misc_Routines.CoFactors.d2r) + \ | |||
|
605 | 2.5*numpy.cos((315.9433 + 2.0*mv - 3.0*me )*Misc_Routines.CoFactors.d2r) + \ | |||
|
606 | 1.6*numpy.cos((345.2533 + 3.0*mv - 4.0*me )*Misc_Routines.CoFactors.d2r) + \ | |||
|
607 | 1.0*numpy.cos((318.15 + 3.0*mv - 5.0*me )*Misc_Routines.CoFactors.d2r) | |||
|
608 | l = l + vencorr | |||
|
609 | ||||
|
610 | # Allow for the Mars perturbations using the mean anomaly of Mars MM | |||
|
611 | mm = 319.529425 + ((19139.858500*t) % 360.0) | |||
|
612 | marscorr = 2.0*numpy.cos((343.8883 - 2.0*mm + 2.0*me)*Misc_Routines.CoFactors.d2r ) + \ | |||
|
613 | 1.8*numpy.cos((200.4017 - 2.0*mm + me)*Misc_Routines.CoFactors.d2r) | |||
|
614 | l = l + marscorr | |||
|
615 | ||||
|
616 | # Allow for the Jupiter perturbations using the mean anomaly of Jupiter MJ | |||
|
617 | mj = 225.328328 + ((3034.6920239*t) % 360.0) | |||
|
618 | jupcorr = 7.2*numpy.cos((179.5317 - mj + me )*Misc_Routines.CoFactors.d2r) + \ | |||
|
619 | 2.6*numpy.cos((263.2167 - mj)*Misc_Routines.CoFactors.d2r) + \ | |||
|
620 | 2.7*numpy.cos((87.1450 - 2.0*mj + 2.0*me)*Misc_Routines.CoFactors.d2r) + \ | |||
|
621 | 1.6*numpy.cos((109.4933 - 2.0*mj + me)*Misc_Routines.CoFactors.d2r) | |||
|
622 | l = l + jupcorr | |||
|
623 | ||||
|
624 | # Allow for Moons perturbations using mean elongation of the Moon from the Sun D | |||
|
625 | d = 350.7376814 + ((445267.11422*t) % 360.0) | |||
|
626 | mooncorr = 6.5*numpy.sin(d*Misc_Routines.CoFactors.d2r) | |||
|
627 | l = l + mooncorr | |||
|
628 | ||||
|
629 | # Allow for long period terms | |||
|
630 | longterm = + 6.4*numpy.sin((231.19 + 20.20*t)*Misc_Routines.CoFactors.d2r) | |||
|
631 | l = l + longterm | |||
|
632 | l = (l + 2592000.0) % 1296000.0 | |||
|
633 | longmed = l/3600.0 | |||
|
634 | ||||
|
635 | # Allow for Aberration | |||
|
636 | l = l - 20.5 | |||
|
637 | ||||
|
638 | # Allow for Nutation using the longitude of the Moons mean node OMEGA | |||
|
639 | omega = 259.183275 - ((1934.142008*t) % 360.0) | |||
|
640 | l = l - 17.2*numpy.sin(omega*Misc_Routines.CoFactors.d2r) | |||
|
641 | ||||
|
642 | # Form the True Obliquity | |||
|
643 | oblt = 23.452294 - 0.0130125*t + (9.2*numpy.cos(omega*Misc_Routines.CoFactors.d2r))/3600.0 | |||
|
644 | ||||
|
645 | # Form Right Ascension and Declination | |||
|
646 | l = l/3600.0 | |||
|
647 | ra = numpy.arctan2((numpy.sin(l*Misc_Routines.CoFactors.d2r)*numpy.cos(oblt*Misc_Routines.CoFactors.d2r)),numpy.cos(l*Misc_Routines.CoFactors.d2r)) | |||
|
648 | ||||
|
649 | neg = numpy.where(ra < 0.0) | |||
|
650 | if neg[0].size > 0: ra[neg] = ra[neg] + 2.0*numpy.pi | |||
|
651 | ||||
|
652 | dec = numpy.arcsin(numpy.sin(l*Misc_Routines.CoFactors.d2r)*numpy.sin(oblt*Misc_Routines.CoFactors.d2r)) | |||
|
653 | ||||
|
654 | if rad==1: | |||
|
655 | oblt = oblt*Misc_Routines.CoFactors.d2r | |||
|
656 | longmed = longmed*Misc_Routines.CoFactors.d2r | |||
|
657 | else: | |||
|
658 | ra = ra/Misc_Routines.CoFactors.d2r | |||
|
659 | dec = dec/Misc_Routines.CoFactors.d2r | |||
|
660 | ||||
|
661 | return ra, dec, longmed, oblt | |||
|
662 | ||||
|
663 | def moonpos(self,jd,rad=0): | |||
|
664 | """ | |||
|
665 | moonpos method computes the RA and Dec of the Moon at specified Julian date(s). | |||
|
666 | ||||
|
667 | Parameters | |||
|
668 | ---------- | |||
|
669 | jd = The julian date of the day (and time), scalar or vector. | |||
|
670 | rad = If this keyword is set and non-zero, then the input and output RAD and DEC | |||
|
671 | vectors are in radian rather than degree. | |||
|
672 | ||||
|
673 | Return | |||
|
674 | ------ | |||
|
675 | ra = The right ascension of the sun at that date in degrees. | |||
|
676 | dec = The declination of the sun at that date in degrees. | |||
|
677 | dist = The Earth-moon distance in kilometers (between the center of the Earth and | |||
|
678 | the center of the moon). | |||
|
679 | geolon = Apparent longitude of the moon in degrees, referred to the ecliptic of the | |||
|
680 | specified date(s). | |||
|
681 | geolat = Apparent latitude the moon in degrees, referred to the ecliptic of the | |||
|
682 | specified date(s). | |||
|
683 | ||||
|
684 | Examples | |||
|
685 | -------- | |||
|
686 | >> jd = 2448724.5 | |||
|
687 | >> [ra,dec,dist,geolon,geolat] = sunpos(jd) | |||
|
688 | >> print ra, dec, dist, geolon, geolat | |||
|
689 | [ 134.68846855] [ 13.76836663] [ 368409.68481613] [ 133.16726428] [-3.22912642] | |||
|
690 | ||||
|
691 | >> [ra,dec,dist,geolon, geolat] = sunpos(jd,rad=1) | |||
|
692 | >> print ra, dec, dist, geolon, geolat | |||
|
693 | [ 2.35075724] [ 0.24030333] [ 368409.68481613] [ 2.32420722] [-0.05635889] | |||
|
694 | ||||
|
695 | >> jd = 2450449.5 + numpy.arange(365) | |||
|
696 | >> [ra,dec,dist,geolon, geolat] = sunpos(jd) | |||
|
697 | ||||
|
698 | Modification history | |||
|
699 | -------------------- | |||
|
700 | Written by Micheal R. Greason, STX Corporation, 31 October 1988. | |||
|
701 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
702 | """ | |||
|
703 | ||||
|
704 | jd = numpy.atleast_1d(jd) | |||
|
705 | ||||
|
706 | # Form time in Julian centuries from 1900. | |||
|
707 | t = (jd - 2451545.0)/36525.0 | |||
|
708 | ||||
|
709 | d_lng = numpy.array([0,2,2,0,0,0,2,2,2,2,0,1,0,2,0,0,4,0,4,2,2,1,1,2,2,4,2,0,2,2,1,2,\ | |||
|
710 | 0,0,2,2,2,4,0,3,2,4,0,2,2,2,4,0,4,1,2,0,1,3,4,2,0,1,2,2]) | |||
|
711 | ||||
|
712 | m_lng = numpy.array([0,0,0,0,1,0,0,-1,0,-1,1,0,1,0,0,0,0,0,0,1,1,0,1,-1,0,0,0,1,0,-1,\ | |||
|
713 | 0,-2,1,2,-2,0,0,-1,0,0,1,-1,2,2,1,-1,0,0,-1,0,1,0,1,0,0,-1,2,1,0,0]) | |||
|
714 | ||||
|
715 | mp_lng = numpy.array([1,-1,0,2,0,0,-2,-1,1,0,-1,0,1,0,1,1,-1,3,-2,-1,0,-1,0,1,2,0,-3,\ | |||
|
716 | -2,-1,-2,1,0,2,0,-1,1,0,-1,2,-1,1,-2,-1,-1,-2,0,1,4,0,-2,0,2,1,-2,-3,2,1,-1,3,-1]) | |||
|
717 | ||||
|
718 | f_lng = numpy.array([0,0,0,0,0,2,0,0,0,0,0,0,0,-2,2,-2,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,\ | |||
|
719 | 0,0,0,0,-2,2,0,2,0,0,0,0,0,0,-2,0,0,0,0,-2,-2,0,0,0,0,0,0,0,-2]) | |||
|
720 | ||||
|
721 | sin_lng = numpy.array([6288774,1274027,658314,213618,-185116,-114332,58793,57066,\ | |||
|
722 | 53322,45758,-40923,-34720,-30383,15327,-12528,10980,10675,10034,8548,-7888,\ | |||
|
723 | -6766,-5163,4987,4036,3994,3861,3665,-2689,-2602,2390,-2348,2236,-2120,-2069,\ | |||
|
724 | 2048,-1773,-1595,1215,-1110,-892,-810,759,-713,-700,691,596,549,537,520,-487,\ | |||
|
725 | -399,-381,351,-340,330,327,-323,299,294,0.0]) | |||
|
726 | ||||
|
727 | cos_lng = numpy.array([-20905355,-3699111,-2955968,-569925,48888,-3149,246158,-152138,\ | |||
|
728 | -170733,-204586,-129620,108743,104755,10321,0,79661,-34782,-23210,-21636,24208,\ | |||
|
729 | 30824,-8379,-16675,-12831,-10445,-11650,14403,-7003,0,10056,6322, -9884,5751,0,\ | |||
|
730 | -4950,4130,0,-3958,0,3258,2616,-1897,-2117,2354,0,0,-1423,-1117,-1571,-1739,0, \ | |||
|
731 | -4421,0,0,0,0,1165,0,0,8752.0]) | |||
|
732 | ||||
|
733 | d_lat = numpy.array([0,0,0,2,2,2,2,0,2,0,2,2,2,2,2,2,2,0,4,0,0,0,1,0,0,0,1,0,4,4,0,4,\ | |||
|
734 | 2,2,2,2,0,2,2,2,2,4,2,2,0,2,1,1,0,2,1,2,0,4,4,1,4,1,4,2]) | |||
|
735 | ||||
|
736 | m_lat = numpy.array([0,0,0,0,0,0,0,0,0,0,-1,0,0,1,-1,-1,-1,1,0,1,0,1,0,1,1,1,0,0,0,0,\ | |||
|
737 | 0,0,0,0,-1,0,0,0,0,1,1,0,-1,-2,0,1,1,1,1,1,0,-1,1,0,-1,0,0,0,-1,-2]) | |||
|
738 | ||||
|
739 | mp_lat = numpy.array([0,1,1,0,-1,-1,0,2,1,2,0,-2,1,0,-1,0,-1,-1,-1,0,0,-1,0,1,1,0,0,\ | |||
|
740 | 3,0,-1,1,-2,0,2,1,-2,3,2,-3,-1,0,0,1,0,1,1,0,0,-2,-1,1,-2,2,-2,-1,1,1,-1,0,0]) | |||
|
741 | ||||
|
742 | f_lat = numpy.array([1,1,-1,-1,1,-1,1,1,-1,-1,-1,-1,1,-1,1,1,-1,-1,-1,1,3,1,1,1,-1,\ | |||
|
743 | -1,-1,1,-1,1,-3,1,-3,-1,-1,1,-1,1,-1,1,1,1,1,-1,3,-1,-1,1,-1,-1,1,-1,1,-1,-1, \ | |||
|
744 | -1,-1,-1,-1,1]) | |||
|
745 | ||||
|
746 | sin_lat = numpy.array([5128122,280602,277693,173237,55413,46271, 32573, 17198, 9266, \ | |||
|
747 | 8822,8216,4324,4200,-3359,2463,2211,2065,-1870,1828,-1794, -1749, -1565, -1491, \ | |||
|
748 | -1475,-1410,-1344,-1335,1107,1021,833,777,671,607,596,491,-451,439,422,421,-366,\ | |||
|
749 | -351,331,315,302,-283,-229,223,223,-220,-220,-185,181,-177,176, 166, -164, 132, \ | |||
|
750 | -119,115,107.0]) | |||
|
751 | ||||
|
752 | # Mean longitude of the moon refered to mean equinox of the date. | |||
|
753 | coeff0 = numpy.array([-1./6.5194e7,1./538841.,-0.0015786,481267.88123421,218.3164477]) | |||
|
754 | lprimed = numpy.poly1d(coeff0) | |||
|
755 | lprimed = lprimed(t) | |||
|
756 | lprimed = self.cirrange(lprimed,rad=0) | |||
|
757 | lprime = lprimed*Misc_Routines.CoFactors.d2r | |||
|
758 | ||||
|
759 | # Mean elongation of the moon | |||
|
760 | coeff1 = numpy.array([-1./1.13065e8,1./545868.,-0.0018819,445267.1114034,297.8501921]) | |||
|
761 | d = numpy.poly1d(coeff1) | |||
|
762 | d = d(t)*Misc_Routines.CoFactors.d2r | |||
|
763 | d = self.cirrange(d,rad=1) | |||
|
764 | ||||
|
765 | # Sun's mean anomaly | |||
|
766 | coeff2 = numpy.array([1.0/2.449e7,-0.0001536,35999.0502909,357.5291092]) | |||
|
767 | M = numpy.poly1d(coeff2) | |||
|
768 | M = M(t)*Misc_Routines.CoFactors.d2r | |||
|
769 | M = self.cirrange(M,rad=1) | |||
|
770 | ||||
|
771 | # Moon's mean anomaly | |||
|
772 | coeff3 = numpy.array([-1.0/1.4712e7,1.0/6.9699e4,0.0087414,477198.8675055,134.9633964]) | |||
|
773 | Mprime = numpy.poly1d(coeff3) | |||
|
774 | Mprime = Mprime(t)*Misc_Routines.CoFactors.d2r | |||
|
775 | Mprime = self.cirrange(Mprime,rad=1) | |||
|
776 | ||||
|
777 | # Moon's argument of latitude | |||
|
778 | coeff4 = numpy.array([1.0/8.6331e8,-1.0/3.526e7,-0.0036539,483202.0175233,93.2720950]) | |||
|
779 | F = numpy.poly1d(coeff4) | |||
|
780 | F = F(t)*Misc_Routines.CoFactors.d2r | |||
|
781 | F = self.cirrange(F,rad=1) | |||
|
782 | ||||
|
783 | # Eccentricity of Earth's orbit around the sun | |||
|
784 | e = 1 - 0.002516*t - 7.4e-6*(t**2.) | |||
|
785 | e2 = e**2. | |||
|
786 | ||||
|
787 | ecorr1 = numpy.where((numpy.abs(m_lng))==1) | |||
|
788 | ecorr2 = numpy.where((numpy.abs(m_lat))==1) | |||
|
789 | ecorr3 = numpy.where((numpy.abs(m_lng))==2) | |||
|
790 | ecorr4 = numpy.where((numpy.abs(m_lat))==2) | |||
|
791 | ||||
|
792 | # Additional arguments. | |||
|
793 | A1 = (119.75 + 131.849*t)*Misc_Routines.CoFactors.d2r | |||
|
794 | A2 = (53.09 + 479264.290*t)*Misc_Routines.CoFactors.d2r | |||
|
795 | A3 = (313.45 + 481266.484*t)*Misc_Routines.CoFactors.d2r | |||
|
796 | suml_add = 3958.*numpy.sin(A1) + 1962.*numpy.sin(lprime - F) + 318*numpy.sin(A2) | |||
|
797 | sumb_add = -2235.*numpy.sin(lprime) + 382.*numpy.sin(A3) + 175.*numpy.sin(A1-F) + \ | |||
|
798 | 175.*numpy.sin(A1 + F) + 127.*numpy.sin(lprime - Mprime) - 115.*numpy.sin(lprime + Mprime) | |||
|
799 | ||||
|
800 | # Sum the periodic terms | |||
|
801 | geolon = numpy.zeros(jd.size) | |||
|
802 | geolat = numpy.zeros(jd.size) | |||
|
803 | dist = numpy.zeros(jd.size) | |||
|
804 | ||||
|
805 | for i in numpy.arange(jd.size): | |||
|
806 | sinlng = sin_lng | |||
|
807 | coslng = cos_lng | |||
|
808 | sinlat = sin_lat | |||
|
809 | ||||
|
810 | sinlng[ecorr1] = e[i]*sinlng[ecorr1] | |||
|
811 | coslng[ecorr1] = e[i]*coslng[ecorr1] | |||
|
812 | sinlat[ecorr2] = e[i]*sinlat[ecorr2] | |||
|
813 | sinlng[ecorr3] = e2[i]*sinlng[ecorr3] | |||
|
814 | coslng[ecorr3] = e2[i]*coslng[ecorr3] | |||
|
815 | sinlat[ecorr4] = e2[i]*sinlat[ecorr4] | |||
|
816 | ||||
|
817 | arg = d_lng*d[i] + m_lng*M[i] + mp_lng*Mprime[i] + f_lng*F[i] | |||
|
818 | geolon[i] = lprimed[i] + (numpy.sum(sinlng*numpy.sin(arg)) + suml_add[i] )/1.e6 | |||
|
819 | dist[i] = 385000.56 + numpy.sum(coslng*numpy.cos(arg))/1.e3 | |||
|
820 | arg = d_lat*d[i] + m_lat*M[i] + mp_lat*Mprime[i] + f_lat*F[i] | |||
|
821 | geolat[i] = (numpy.sum(sinlat*numpy.sin(arg)) + sumb_add[i])/1.e6 | |||
|
822 | ||||
|
823 | [nlon, elon] = self.nutate(jd) | |||
|
824 | geolon = geolon + nlon/3.6e3 | |||
|
825 | geolon = self.cirrange(geolon,rad=0) | |||
|
826 | lamb = geolon*Misc_Routines.CoFactors.d2r | |||
|
827 | beta = geolat*Misc_Routines.CoFactors.d2r | |||
|
828 | ||||
|
829 | # Find mean obliquity and convert lamb, beta to RA, Dec | |||
|
830 | c = numpy.array([2.45,5.79,27.87,7.12,-39.05,-249.67,-51.38,1999.25,-1.55,-4680.93, \ | |||
|
831 | 21.448]) | |||
|
832 | junk = numpy.poly1d(c); | |||
|
833 | epsilon = 23. + (26./60.) + (junk(t/1.e2)/3600.) | |||
|
834 | # True obliquity in radians | |||
|
835 | eps = (epsilon + elon/3600. )*Misc_Routines.CoFactors.d2r | |||
|
836 | ||||
|
837 | ra = numpy.arctan2(numpy.sin(lamb)*numpy.cos(eps)-numpy.tan(beta)*numpy.sin(eps),numpy.cos(lamb)) | |||
|
838 | ra = self.cirrange(ra,rad=1) | |||
|
839 | ||||
|
840 | dec = numpy.arcsin(numpy.sin(beta)*numpy.cos(eps) + numpy.cos(beta)*numpy.sin(eps)*numpy.sin(lamb)) | |||
|
841 | ||||
|
842 | if rad==1: | |||
|
843 | geolon = lamb | |||
|
844 | geolat = beta | |||
|
845 | else: | |||
|
846 | ra = ra/Misc_Routines.CoFactors.d2r | |||
|
847 | dec = dec/Misc_Routines.CoFactors.d2r | |||
|
848 | ||||
|
849 | return ra, dec, dist, geolon, geolat | |||
|
850 | ||||
|
851 | def hydrapos(self): | |||
|
852 | """ | |||
|
853 | hydrapos method returns RA and Dec provided by Bill Coles (Oct 2003). | |||
|
854 | ||||
|
855 | Parameters | |||
|
856 | ---------- | |||
|
857 | None | |||
|
858 | ||||
|
859 | Return | |||
|
860 | ------ | |||
|
861 | ra = The right ascension of the sun at that date in degrees. | |||
|
862 | dec = The declination of the sun at that date in degrees. | |||
|
863 | Examples | |||
|
864 | -------- | |||
|
865 | >> [ra,dec] = hydrapos() | |||
|
866 | >> print ra, dec | |||
|
867 | 139.45 -12.0833333333 | |||
|
868 | ||||
|
869 | Modification history | |||
|
870 | -------------------- | |||
|
871 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
872 | """ | |||
|
873 | ||||
|
874 | ra = (9. + 17.8/60.)*15. | |||
|
875 | dec = -(12. + 5./60.) | |||
|
876 | ||||
|
877 | return ra, dec | |||
|
878 | ||||
|
879 | ||||
|
880 | def skynoise_jro(self,dec_cut=-11.95,filename='skynoise_jro.dat',filepath=None): | |||
|
881 | """ | |||
|
882 | hydrapos returns RA and Dec provided by Bill Coles (Oct 2003). | |||
|
883 | ||||
|
884 | Parameters | |||
|
885 | ---------- | |||
|
886 | dec_cut = A scalar giving the declination to get a cut of the skynoise over Jica- | |||
|
887 | marca. The default value is -11.95. | |||
|
888 | filename = A string to specify name the skynoise file. The default value is skynoi- | |||
|
889 | se_jro.dat | |||
|
890 | ||||
|
891 | Return | |||
|
892 | ------ | |||
|
893 | maxra = The maximum right ascension to the declination used to get a cut. | |||
|
894 | ra = The right ascension. | |||
|
895 | Examples | |||
|
896 | -------- | |||
|
897 | >> [maxra,ra] = skynoise_jro() | |||
|
898 | >> print maxra, ra | |||
|
899 | 139.45 -12.0833333333 | |||
|
900 | ||||
|
901 | Modification history | |||
|
902 | -------------------- | |||
|
903 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
904 | """ | |||
|
905 | ||||
|
906 | if filepath==None: | |||
|
907 | filepath = '/app/utils/' | |||
|
908 | ||||
|
909 | f = open(os.path.join(filepath,filename),'rb') | |||
|
910 | ||||
|
911 | # Reading SkyNoise Power (lineal scale) | |||
|
912 | ha_sky = numpy.fromfile(f,numpy.dtype([('var','<f4')]),480*20) | |||
|
913 | ha_sky = ha_sky['var'].reshape(20,480).transpose() | |||
|
914 | ||||
|
915 | dec_sky = numpy.fromfile(f,numpy.dtype([('var','<f4')]),480*20) | |||
|
916 | dec_sky = dec_sky['var'].reshape((20,480)).transpose() | |||
|
917 | ||||
|
918 | tmp_sky = numpy.fromfile(f,numpy.dtype([('var','<f4')]),480*20) | |||
|
919 | tmp_sky = tmp_sky['var'].reshape((20,480)).transpose() | |||
|
920 | ||||
|
921 | f.close() | |||
|
922 | ||||
|
923 | nha = 480 | |||
|
924 | tmp_cut = numpy.zeros(nha) | |||
|
925 | for iha in numpy.arange(nha): | |||
|
926 | tck = scipy.interpolate.splrep(dec_sky[iha,:],tmp_sky[iha,:],s=0) | |||
|
927 | tmp_cut[iha] = scipy.interpolate.splev(dec_cut,tck,der=0) | |||
|
928 | ||||
|
929 | ptr = numpy.nanargmax(tmp_cut) | |||
|
930 | ||||
|
931 | maxra = ha_sky[ptr,0] | |||
|
932 | ra = ha_sky[:,0] | |||
|
933 | ||||
|
934 | return maxra, ra | |||
|
935 | ||||
|
936 | def skyNoise(self,jd,ut=-5.0,longitude=-76.87,filename='galaxy.txt',filepath=None): | |||
|
937 | """ | |||
|
938 | hydrapos returns RA and Dec provided by Bill Coles (Oct 2003). | |||
|
939 | ||||
|
940 | Parameters | |||
|
941 | ---------- | |||
|
942 | jd = The julian date of the day (and time), scalar or vector. | |||
|
943 | ||||
|
944 | dec_cut = A scalar giving the declination to get a cut of the skynoise over Jica- | |||
|
945 | marca. The default value is -11.95. | |||
|
946 | filename = A string to specify name the skynoise file. The default value is skynoi- | |||
|
947 | se_jro.dat | |||
|
948 | ||||
|
949 | Return | |||
|
950 | ------ | |||
|
951 | maxra = The maximum right ascension to the declination used to get a cut. | |||
|
952 | ra = The right ascension. | |||
|
953 | ||||
|
954 | Examples | |||
|
955 | -------- | |||
|
956 | >> [maxra,ra] = skynoise_jro() | |||
|
957 | >> print maxra, ra | |||
|
958 | 139.45 -12.0833333333 | |||
|
959 | ||||
|
960 | Modification history | |||
|
961 | -------------------- | |||
|
962 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
963 | """ | |||
|
964 | ||||
|
965 | # Defining date to compute SkyNoise. | |||
|
966 | [year, month, dom, hour, mis, secs] = TimeTools.Julian(jd).change2time() | |||
|
967 | is_dom = (month==9) & (dom==21) | |||
|
968 | if is_dom: | |||
|
969 | tmp = jd | |||
|
970 | jd = TimeTools.Time(year,9,22).change2julian() | |||
|
971 | dom = 22 | |||
|
972 | ||||
|
973 | # Reading SkyNoise | |||
|
974 | if filepath==None:filepath='./resource' | |||
|
975 | f = open(os.path.join(filepath,filename)) | |||
|
976 | ||||
|
977 | lines = f.read() | |||
|
978 | f.close() | |||
|
979 | ||||
|
980 | nlines = 99 | |||
|
981 | lines = lines.split('\n') | |||
|
982 | data = numpy.zeros((2,nlines))*numpy.float32(0.) | |||
|
983 | for ii in numpy.arange(nlines): | |||
|
984 | line = numpy.array([lines[ii][0:6],lines[ii][6:]]) | |||
|
985 | data[:,ii] = numpy.float32(line) | |||
|
986 | ||||
|
987 | # Getting SkyNoise to the date desired. | |||
|
988 | otime = data[0,:]*60.0 | |||
|
989 | opowr = data[1,:] | |||
|
990 | ||||
|
991 | hour = numpy.array([0,23]); | |||
|
992 | mins = numpy.array([0,59]); | |||
|
993 | secs = numpy.array([0,59]); | |||
|
994 | LTrange = TimeTools.Time(year,month,dom,hour,mins,secs).change2julday() | |||
|
995 | LTtime = LTrange[0] + numpy.arange(1440)*((LTrange[1] - LTrange[0])/(1440.-1)) | |||
|
996 | lst = TimeTools.Julian(LTtime + (-3600.*ut/86400.)).change2lst() | |||
|
997 | ||||
|
998 | ipowr = lst*0.0 | |||
|
999 | # Interpolating using scipy (inside max and min "x") | |||
|
1000 | otime = otime/3600. | |||
|
1001 | val = numpy.where((lst>numpy.min(otime)) & (lst<numpy.max(otime))); val = val[0] | |||
|
1002 | tck = scipy.interpolate.interp1d(otime,opowr) | |||
|
1003 | ipowr[val] = tck(lst[val]) | |||
|
1004 | ||||
|
1005 | # Extrapolating above maximum time data (23.75). | |||
|
1006 | uval = numpy.where(lst>numpy.max(otime)) | |||
|
1007 | if uval[0].size>0: | |||
|
1008 | ii = numpy.min(uval[0]) | |||
|
1009 | m = (ipowr[ii-1] - ipowr[ii-2])/(lst[ii-1] - lst[ii-2]) | |||
|
1010 | b = ipowr[ii-1] - m*lst[ii-1] | |||
|
1011 | ipowr[uval] = m*lst[uval] + b | |||
|
1012 | ||||
|
1013 | if is_dom: | |||
|
1014 | lst = numpy.roll(lst,4) | |||
|
1015 | ipowr = numpy.roll(ipowr,4) | |||
|
1016 | ||||
|
1017 | new_lst = numpy.int32(lst*3600.) | |||
|
1018 | new_pow = ipowr | |||
|
1019 | ||||
|
1020 | return ipowr, LTtime, lst | |||
|
1021 | ||||
|
1022 | ||||
|
1023 | class AltAz(EquatorialCorrections): | |||
|
1024 | def __init__(self,alt,az,jd,lat=-11.95,lon=-76.8667,WS=0,altitude=500,nutate_=0,precess_=0,\ | |||
|
1025 | aberration_=0,B1950=0): | |||
|
1026 | """ | |||
|
1027 | The AltAz class creates an object which represents the target position in horizontal | |||
|
1028 | coordinates (alt-az) and allows to convert (using the methods) from this coordinate | |||
|
1029 | system to others (e.g. Equatorial). | |||
|
1030 | ||||
|
1031 | Parameters | |||
|
1032 | ---------- | |||
|
1033 | alt = Altitude in degrees. Scalar or vector. | |||
|
1034 | az = Azimuth angle in degrees (measured EAST from NORTH, but see keyword WS). Sca- | |||
|
1035 | lar or vector. | |||
|
1036 | jd = Julian date. Scalar or vector. | |||
|
1037 | lat = North geodetic latitude of location in degrees. The default value is -11.95. | |||
|
1038 | lon = East longitude of location in degrees. The default value is -76.8667. | |||
|
1039 | WS = Set this to 1 to get the azimuth measured westward from south. | |||
|
1040 | altitude = The altitude of the observing location, in meters. The default 500. | |||
|
1041 | nutate_ = Set this to 1 to force nutation, 0 for no nutation. | |||
|
1042 | precess_ = Set this to 1 to force precession, 0 for no precession. | |||
|
1043 | aberration_ = Set this to 1 to force aberration correction, 0 for no correction. | |||
|
1044 | B1950 = Set this if your RA and DEC are specified in B1950, FK4 coordinates (ins- | |||
|
1045 | tead of J2000, FK5) | |||
|
1046 | ||||
|
1047 | Modification History | |||
|
1048 | -------------------- | |||
|
1049 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 26 September 2009. | |||
|
1050 | """ | |||
|
1051 | ||||
|
1052 | EquatorialCorrections.__init__(self) | |||
|
1053 | ||||
|
1054 | self.alt = numpy.atleast_1d(alt) | |||
|
1055 | self.az = numpy.atleast_1d(az) | |||
|
1056 | self.jd = numpy.atleast_1d(jd) | |||
|
1057 | self.lat = lat | |||
|
1058 | self.lon = lon | |||
|
1059 | self.WS = WS | |||
|
1060 | self.altitude = altitude | |||
|
1061 | ||||
|
1062 | self.nutate_ = nutate_ | |||
|
1063 | self.aberration_ = aberration_ | |||
|
1064 | self.precess_ = precess_ | |||
|
1065 | self.B1950 = B1950 | |||
|
1066 | ||||
|
1067 | def change2equatorial(self): | |||
|
1068 | """ | |||
|
1069 | change2equatorial method converts horizon (Alt-Az) coordinates to equatorial coordi- | |||
|
1070 | nates (ra-dec). | |||
|
1071 | ||||
|
1072 | Return | |||
|
1073 | ------ | |||
|
1074 | ra = Right ascension of object (J2000) in degrees (FK5). Scalar or vector. | |||
|
1075 | dec = Declination of object (J2000), in degrees (FK5). Scalar or vector. | |||
|
1076 | ha = Hour angle in degrees. | |||
|
1077 | ||||
|
1078 | Example | |||
|
1079 | ------- | |||
|
1080 | >> alt = 88.5401 | |||
|
1081 | >> az = -128.990 | |||
|
1082 | >> jd = 2452640.5 | |||
|
1083 | >> ObjAltAz = AltAz(alt,az,jd) | |||
|
1084 | >> [ra, dec, ha] = ObjAltAz.change2equatorial() | |||
|
1085 | >> print ra, dec, ha | |||
|
1086 | [ 22.20280632] [-12.86610025] [ 1.1638927] | |||
|
1087 | ||||
|
1088 | Modification History | |||
|
1089 | -------------------- | |||
|
1090 | Written Chris O'Dell Univ. of Wisconsin-Madison, May 2002. | |||
|
1091 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
1092 | """ | |||
|
1093 | ||||
|
1094 | az = self.az | |||
|
1095 | alt = self.alt | |||
|
1096 | if self.WS>0:az = az -180. | |||
|
1097 | ra_tmp = numpy.zeros(numpy.size(self.jd)) + 45. | |||
|
1098 | dec_tmp = numpy.zeros(numpy.size(self.jd)) + 45. | |||
|
1099 | [dra1,ddec1,eps,d_psi,d_eps] = self.co_nutate(self.jd,ra_tmp, dec_tmp) | |||
|
1100 | ||||
|
1101 | # Getting local mean sidereal time (lmst) | |||
|
1102 | lmst = TimeTools.Julian(self.jd[0]).change2lst() | |||
|
1103 | lmst = lmst*Misc_Routines.CoFactors.h2d | |||
|
1104 | # Getting local apparent sidereal time (last) | |||
|
1105 | last = lmst + d_psi*numpy.cos(eps)/3600. | |||
|
1106 | ||||
|
1107 | # Now do the spherical trig to get APPARENT hour angle and declination (Degrees). | |||
|
1108 | [ha, dec] = self.change2HaDec() | |||
|
1109 | ||||
|
1110 | # Finding Right Ascension (in degrees, from 0 to 360.) | |||
|
1111 | ra = (last - ha + 360.) % 360. | |||
|
1112 | ||||
|
1113 | # Calculate NUTATION and ABERRATION Correction to Ra-Dec | |||
|
1114 | [dra1, ddec1,eps,d_psi,d_eps] = self.co_nutate(self.jd,ra,dec) | |||
|
1115 | [dra2,ddec2,eps] = self.co_aberration(self.jd,ra,dec) | |||
|
1116 | ||||
|
1117 | # Make Nutation and Aberration correction (if wanted) | |||
|
1118 | ra = ra - (dra1*self.nutate_ + dra2*self.aberration_)/3600. | |||
|
1119 | dec = dec - (ddec1*self.nutate_ + ddec2*self.aberration_)/3600. | |||
|
1120 | ||||
|
1121 | # Computing current equinox | |||
|
1122 | j_now = (self.jd - 2451545.)/365.25 + 2000 | |||
|
1123 | ||||
|
1124 | # Precess coordinates to current date | |||
|
1125 | if self.precess_==1: | |||
|
1126 | njd = numpy.size(self.jd) | |||
|
1127 | for ii in numpy.arange(njd): | |||
|
1128 | ra_i = ra[ii] | |||
|
1129 | dec_i = dec[ii] | |||
|
1130 | now = j_now[ii] | |||
|
1131 | ||||
|
1132 | if self.B1950==1: | |||
|
1133 | [ra_i,dec_i] = self.precess(ra_i,dec_i,now,1950.,FK4=1) | |||
|
1134 | elif self.B1950==0: | |||
|
1135 | [ra_i,dec_i] = self.precess(ra_i,dec_i,now,2000.,FK4=0) | |||
|
1136 | ||||
|
1137 | ra[ii] = ra_i | |||
|
1138 | dec[ii] = dec_i | |||
|
1139 | ||||
|
1140 | return ra, dec, ha | |||
|
1141 | ||||
|
1142 | def change2HaDec(self): | |||
|
1143 | """ | |||
|
1144 | change2HaDec method converts from horizon (Alt-Az) coordinates to hour angle and de- | |||
|
1145 | clination. | |||
|
1146 | ||||
|
1147 | Return | |||
|
1148 | ------ | |||
|
1149 | ha = The local apparent hour angle, in degrees. The hour angle is the time that ri- | |||
|
1150 | ght ascension of 0 hours crosses the local meridian. It is unambiguisoly defined. | |||
|
1151 | dec = The local apparent declination, in degrees. | |||
|
1152 | ||||
|
1153 | Example | |||
|
1154 | ------- | |||
|
1155 | >> alt = 88.5401 | |||
|
1156 | >> az = -128.990 | |||
|
1157 | >> jd = 2452640.5 | |||
|
1158 | >> ObjAltAz = AltAz(alt,az,jd) | |||
|
1159 | >> [ha, dec] = ObjAltAz.change2HaDec() | |||
|
1160 | >> print ha, dec | |||
|
1161 | [ 1.1638927] [-12.86610025] | |||
|
1162 | ||||
|
1163 | Modification History | |||
|
1164 | -------------------- | |||
|
1165 | Written Chris O'Dell Univ. of Wisconsin-Madison, May 2002. | |||
|
1166 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
1167 | """ | |||
|
1168 | ||||
|
1169 | alt_r = numpy.atleast_1d(self.alt*Misc_Routines.CoFactors.d2r) | |||
|
1170 | az_r = numpy.atleast_1d(self.az*Misc_Routines.CoFactors.d2r) | |||
|
1171 | lat_r = numpy.atleast_1d(self.lat*Misc_Routines.CoFactors.d2r) | |||
|
1172 | ||||
|
1173 | # Find local hour angle (in degrees, from 0 to 360.) | |||
|
1174 | y_ha = -1*numpy.sin(az_r)*numpy.cos(alt_r) | |||
|
1175 | x_ha = -1*numpy.cos(az_r)*numpy.sin(lat_r)*numpy.cos(alt_r) + numpy.sin(alt_r)*numpy.cos(lat_r) | |||
|
1176 | ||||
|
1177 | ha = numpy.arctan2(y_ha,x_ha) | |||
|
1178 | ha = ha/Misc_Routines.CoFactors.d2r | |||
|
1179 | ||||
|
1180 | w = numpy.where(ha<0.) | |||
|
1181 | if w[0].size>0:ha[w] = ha[w] + 360. | |||
|
1182 | ha = ha % 360. | |||
|
1183 | ||||
|
1184 | # Find declination (positive if north of celestial equatorial, negative if south) | |||
|
1185 | sindec = numpy.sin(lat_r)*numpy.sin(alt_r) + numpy.cos(lat_r)*numpy.cos(alt_r)*numpy.cos(az_r) | |||
|
1186 | dec = numpy.arcsin(sindec)/Misc_Routines.CoFactors.d2r | |||
|
1187 | ||||
|
1188 | return ha, dec | |||
|
1189 | ||||
|
1190 | ||||
|
1191 | class Equatorial(EquatorialCorrections): | |||
|
1192 | def __init__(self,ra,dec,jd,lat=-11.95,lon=-76.8667,WS=0,altitude=500,nutate_=0,precess_=0,\ | |||
|
1193 | aberration_=0,B1950=0): | |||
|
1194 | """ | |||
|
1195 | The Equatorial class creates an object which represents the target position in equa- | |||
|
1196 | torial coordinates (ha-dec) and allows to convert (using the class methods) from | |||
|
1197 | this coordinate system to others (e.g. AltAz). | |||
|
1198 | ||||
|
1199 | Parameters | |||
|
1200 | ---------- | |||
|
1201 | ra = Right ascension of object (J2000) in degrees (FK5). Scalar or vector. | |||
|
1202 | dec = Declination of object (J2000), in degrees (FK5). Scalar or vector. | |||
|
1203 | jd = Julian date. Scalar or vector. | |||
|
1204 | lat = North geodetic latitude of location in degrees. The default value is -11.95. | |||
|
1205 | lon = East longitude of location in degrees. The default value is -76.8667. | |||
|
1206 | WS = Set this to 1 to get the azimuth measured westward from south. | |||
|
1207 | altitude = The altitude of the observing location, in meters. The default 500. | |||
|
1208 | nutate = Set this to 1 to force nutation, 0 for no nutation. | |||
|
1209 | precess = Set this to 1 to force precession, 0 for no precession. | |||
|
1210 | aberration = Set this to 1 to force aberration correction, 0 for no correction. | |||
|
1211 | B1950 = Set this if your RA and DEC are specified in B1950, FK4 coordinates (ins- | |||
|
1212 | tead of J2000, FK5) | |||
|
1213 | ||||
|
1214 | Modification History | |||
|
1215 | -------------------- | |||
|
1216 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 29 September 2009. | |||
|
1217 | """ | |||
|
1218 | ||||
|
1219 | EquatorialCorrections.__init__(self) | |||
|
1220 | ||||
|
1221 | self.ra = numpy.atleast_1d(ra) | |||
|
1222 | self.dec = numpy.atleast_1d(dec) | |||
|
1223 | self.jd = numpy.atleast_1d(jd) | |||
|
1224 | self.lat = lat | |||
|
1225 | self.lon = lon | |||
|
1226 | self.WS = WS | |||
|
1227 | self.altitude = altitude | |||
|
1228 | ||||
|
1229 | self.nutate_ = nutate_ | |||
|
1230 | self.aberration_ = aberration_ | |||
|
1231 | self.precess_ = precess_ | |||
|
1232 | self.B1950 = B1950 | |||
|
1233 | ||||
|
1234 | def change2AltAz(self): | |||
|
1235 | """ | |||
|
1236 | change2AltAz method converts from equatorial coordinates (ha-dec) to horizon coordi- | |||
|
1237 | nates (alt-az). | |||
|
1238 | ||||
|
1239 | Return | |||
|
1240 | ------ | |||
|
1241 | alt = Altitude in degrees. Scalar or vector. | |||
|
1242 | az = Azimuth angle in degrees (measured EAST from NORTH, but see keyword WS). Sca- | |||
|
1243 | lar or vector. | |||
|
1244 | ha = Hour angle in degrees. | |||
|
1245 | ||||
|
1246 | Example | |||
|
1247 | ------- | |||
|
1248 | >> ra = 43.370609 | |||
|
1249 | >> dec = -28.0000 | |||
|
1250 | >> jd = 2452640.5 | |||
|
1251 | >> ObjEq = Equatorial(ra,dec,jd) | |||
|
1252 | >> [alt, az, ha] = ObjEq.change2AltAz() | |||
|
1253 | >> print alt, az, ha | |||
|
1254 | [ 65.3546497] [ 133.58753124] [ 339.99609002] | |||
|
1255 | ||||
|
1256 | Modification History | |||
|
1257 | -------------------- | |||
|
1258 | Written Chris O'Dell Univ. of Wisconsin-Madison. May 2002 | |||
|
1259 | Converted to Python by Freddy R. Galindo, ROJ, 29 September 2009. | |||
|
1260 | """ | |||
|
1261 | ||||
|
1262 | ra = self.ra | |||
|
1263 | dec = self.dec | |||
|
1264 | ||||
|
1265 | # Computing current equinox | |||
|
1266 | j_now = (self.jd - 2451545.)/365.25 + 2000 | |||
|
1267 | ||||
|
1268 | # Precess coordinates to current date | |||
|
1269 | if self.precess_==1: | |||
|
1270 | njd = numpy.size(self.jd) | |||
|
1271 | for ii in numpy.arange(njd): | |||
|
1272 | ra_i = ra[ii] | |||
|
1273 | dec_i = dec[ii] | |||
|
1274 | now = j_now[ii] | |||
|
1275 | ||||
|
1276 | if self.B1950==1: | |||
|
1277 | [ra_i,dec_i] = self.precess(ra_i,dec_i,now,1950.,FK4=1) | |||
|
1278 | elif self.B1950==0: | |||
|
1279 | [ra_i,dec_i] = self.precess(ra_i,dec_i,now,2000.,FK4=0) | |||
|
1280 | ||||
|
1281 | ra[ii] = ra_i | |||
|
1282 | dec[ii] = dec_i | |||
|
1283 | ||||
|
1284 | # Calculate NUTATION and ABERRATION Correction to Ra-Dec | |||
|
1285 | [dra1, ddec1,eps,d_psi,d_eps] = self.co_nutate(self.jd,ra,dec) | |||
|
1286 | [dra2,ddec2,eps] = self.co_aberration(self.jd,ra,dec) | |||
|
1287 | ||||
|
1288 | # Make Nutation and Aberration correction (if wanted) | |||
|
1289 | ra = ra + (dra1*self.nutate_ + dra2*self.aberration_)/3600. | |||
|
1290 | dec = dec + (ddec1*self.nutate_ + ddec2*self.aberration_)/3600. | |||
|
1291 | ||||
|
1292 | # Getting local mean sidereal time (lmst) | |||
|
1293 | lmst = TimeTools.Julian(self.jd).change2lst() | |||
|
1294 | ||||
|
1295 | lmst = lmst*Misc_Routines.CoFactors.h2d | |||
|
1296 | # Getting local apparent sidereal time (last) | |||
|
1297 | last = lmst + d_psi*numpy.cos(eps)/3600. | |||
|
1298 | ||||
|
1299 | # Finding Hour Angle (in degrees, from 0 to 360.) | |||
|
1300 | ha = last - ra | |||
|
1301 | w = numpy.where(ha<0.) | |||
|
1302 | if w[0].size>0:ha[w] = ha[w] + 360. | |||
|
1303 | ha = ha % 360. | |||
|
1304 | ||||
|
1305 | # Now do the spherical trig to get APPARENT hour angle and declination (Degrees). | |||
|
1306 | [alt, az] = self.HaDec2AltAz(ha,dec) | |||
|
1307 | ||||
|
1308 | return alt, az, ha | |||
|
1309 | ||||
|
1310 | def HaDec2AltAz(self,ha,dec): | |||
|
1311 | """ | |||
|
1312 | HaDec2AltAz convert hour angle and declination (ha-dec) to horizon coords (alt-az). | |||
|
1313 | ||||
|
1314 | Parameters | |||
|
1315 | ---------- | |||
|
1316 | ha = The local apparent hour angle, in DEGREES, scalar or vector. | |||
|
1317 | dec = The local apparent declination, in DEGREES, scalar or vector. | |||
|
1318 | ||||
|
1319 | Return | |||
|
1320 | ------ | |||
|
1321 | alt = Altitude in degrees. Scalar or vector. | |||
|
1322 | az = Azimuth angle in degrees (measured EAST from NORTH, but see keyword WS). Sca- | |||
|
1323 | lar or vector. | |||
|
1324 | ||||
|
1325 | Modification History | |||
|
1326 | -------------------- | |||
|
1327 | Written Chris O'Dell Univ. of Wisconsin-Madison, May 2002. | |||
|
1328 | Converted to Python by Freddy R. Galindo, ROJ, 26 September 2009. | |||
|
1329 | """ | |||
|
1330 | ||||
|
1331 | sh = numpy.sin(ha*Misc_Routines.CoFactors.d2r) ; ch = numpy.cos(ha*Misc_Routines.CoFactors.d2r) | |||
|
1332 | sd = numpy.sin(dec*Misc_Routines.CoFactors.d2r) ; cd = numpy.cos(dec*Misc_Routines.CoFactors.d2r) | |||
|
1333 | sl = numpy.sin(self.lat*Misc_Routines.CoFactors.d2r) ; cl = numpy.cos(self.lat*Misc_Routines.CoFactors.d2r) | |||
|
1334 | ||||
|
1335 | x = -1*ch*cd*sl + sd*cl | |||
|
1336 | y = -1*sh*cd | |||
|
1337 | z = ch*cd*cl + sd*sl | |||
|
1338 | r = numpy.sqrt(x**2. + y**2.) | |||
|
1339 | ||||
|
1340 | az = numpy.arctan2(y,x)/Misc_Routines.CoFactors.d2r | |||
|
1341 | alt = numpy.arctan2(z,r)/Misc_Routines.CoFactors.d2r | |||
|
1342 | ||||
|
1343 | # correct for negative az. | |||
|
1344 | w = numpy.where(az<0.) | |||
|
1345 | if w[0].size>0:az[w] = az[w] + 360. | |||
|
1346 | ||||
|
1347 | # Convert az to West from South, if desired | |||
|
1348 | if self.WS==1: az = (az + 180.) % 360. | |||
|
1349 | ||||
|
1350 | return alt, az | |||
|
1351 | ||||
|
1352 | ||||
|
1353 | class Geodetic(): | |||
|
1354 | def __init__(self,lat=-11.95,alt=0): | |||
|
1355 | """ | |||
|
1356 | The Geodetic class creates an object which represents the real position on earth of | |||
|
1357 | a target (Geodetic Coordinates: lat-alt) and allows to convert (using the class me- | |||
|
1358 | thods) from this coordinate system to others (e.g. geocentric). | |||
|
1359 | ||||
|
1360 | Parameters | |||
|
1361 | ---------- | |||
|
1362 | lat = Geodetic latitude of location in degrees. The default value is -11.95. | |||
|
1363 | ||||
|
1364 | alt = Geodetic altitude (km). The default value is 0. | |||
|
1365 | ||||
|
1366 | Modification History | |||
|
1367 | -------------------- | |||
|
1368 | Converted to Object-oriented Programming by Freddy R. Galindo, ROJ, 02 October 2009. | |||
|
1369 | """ | |||
|
1370 | ||||
|
1371 | self.lat = numpy.atleast_1d(lat) | |||
|
1372 | self.alt = numpy.atleast_1d(alt) | |||
|
1373 | ||||
|
1374 | self.a = 6378.16 | |||
|
1375 | self.ab2 = 1.0067397 | |||
|
1376 | self.ep2 = 0.0067397 | |||
|
1377 | ||||
|
1378 | def change2geocentric(self): | |||
|
1379 | """ | |||
|
1380 | change2geocentric method converts from Geodetic to Geocentric coordinates. The re- | |||
|
1381 | ference geoid is that adopted by the IAU in 1964. | |||
|
1382 | ||||
|
1383 | Return | |||
|
1384 | ------ | |||
|
1385 | gclat = Geocentric latitude (in degrees), scalar or vector. | |||
|
1386 | gcalt = Geocentric radial distance (km), scalar or vector. | |||
|
1387 | ||||
|
1388 | Example | |||
|
1389 | ------- | |||
|
1390 | >> ObjGeoid = Geodetic(lat=-11.95,alt=0) | |||
|
1391 | >> [gclat, gcalt] = ObjGeoid.change2geocentric() | |||
|
1392 | >> print gclat, gcalt | |||
|
1393 | [-11.87227742] [ 6377.25048195] | |||
|
1394 | ||||
|
1395 | Modification History | |||
|
1396 | -------------------- | |||
|
1397 | Converted to Python by Freddy R. Galindo, ROJ, 02 October 2009. | |||
|
1398 | """ | |||
|
1399 | ||||
|
1400 | gdl = self.lat*Misc_Routines.CoFactors.d2r | |||
|
1401 | slat = numpy.sin(gdl) | |||
|
1402 | clat = numpy.cos(gdl) | |||
|
1403 | slat2 = slat**2. | |||
|
1404 | clat2 = (self.ab2*clat)**2. | |||
|
1405 | ||||
|
1406 | sbet = slat/numpy.sqrt(slat2 + clat2) | |||
|
1407 | sbet2 = (sbet**2.) # < 1 | |||
|
1408 | noval = numpy.where(sbet2>1) | |||
|
1409 | if noval[0].size>0:sbet2[noval] = 1 | |||
|
1410 | cbet = numpy.sqrt(1. - sbet2) | |||
|
1411 | ||||
|
1412 | rgeoid = self.a/numpy.sqrt(1. + self.ep2*sbet2) | |||
|
1413 | ||||
|
1414 | x = rgeoid*cbet + self.alt*clat | |||
|
1415 | y = rgeoid*sbet + self.alt*slat | |||
|
1416 | ||||
|
1417 | gcalt = numpy.sqrt(x**2. + y**2.) | |||
|
1418 | gclat = numpy.arctan2(y,x)/Misc_Routines.CoFactors.d2r | |||
|
1419 | ||||
|
1420 | return gclat, gcalt |
@@ -0,0 +1,61 | |||||
|
1 | """ | |||
|
2 | The module MISC_ROUTINES gathers classes and functions which are useful for daily processing. As an | |||
|
3 | example we have conversion factor or universal constants. | |||
|
4 | ||||
|
5 | MODULES CALLED: | |||
|
6 | NUMPY, SYS | |||
|
7 | ||||
|
8 | MODIFICATION HISTORY: | |||
|
9 | Created by Ing. Freddy Galindo (frederickgalindo@gmail.com). ROJ, 21 October 2009. | |||
|
10 | """ | |||
|
11 | ||||
|
12 | import numpy | |||
|
13 | import sys | |||
|
14 | ||||
|
15 | class CoFactors(): | |||
|
16 | """ | |||
|
17 | CoFactor class used to call pre-defined conversion factor (e.g. degree to radian). The cu- | |||
|
18 | The current available factor are: | |||
|
19 | ||||
|
20 | d2r = degree to radian. | |||
|
21 | s2r = seconds to radian?, degree to arcsecond.? | |||
|
22 | h2r = hour to radian. | |||
|
23 | h2d = hour to degree | |||
|
24 | """ | |||
|
25 | ||||
|
26 | d2r = numpy.pi/180. | |||
|
27 | s2r = numpy.pi/(180.*3600.) | |||
|
28 | h2r = numpy.pi/12. | |||
|
29 | h2d = 15. | |||
|
30 | ||||
|
31 | ||||
|
32 | class Vector: | |||
|
33 | """ | |||
|
34 | direction = 0 Polar to rectangular; direction=1 rectangular to polar | |||
|
35 | """ | |||
|
36 | def __init__(self,vect,direction=0): | |||
|
37 | nsize = numpy.size(vect) | |||
|
38 | if nsize <= 3: | |||
|
39 | vect = vect.reshape(1,nsize) | |||
|
40 | ||||
|
41 | self.vect = vect | |||
|
42 | self.dirc = direction | |||
|
43 | ||||
|
44 | ||||
|
45 | ||||
|
46 | def Polar2Rect(self): | |||
|
47 | if self.dirc == 0: | |||
|
48 | jvect = self.vect*numpy.pi/180. | |||
|
49 | mmx = numpy.cos(jvect[:,1])*numpy.sin(jvect[:,0]) | |||
|
50 | mmy = numpy.cos(jvect[:,1])*numpy.cos(jvect[:,0]) | |||
|
51 | mmz = numpy.sin(jvect[:,1]) | |||
|
52 | mm = numpy.array([mmx,mmy,mmz]).transpose() | |||
|
53 | ||||
|
54 | elif self.dirc == 1: | |||
|
55 | mm = [numpy.arctan2(self.vect[:,0],self.vect[:,1]),numpy.arcsin(self.vect[:,2])] | |||
|
56 | mm = numpy.array(mm)*180./numpy.pi | |||
|
57 | ||||
|
58 | return mm | |||
|
59 | ||||
|
60 | ||||
|
61 | No newline at end of file |
@@ -0,0 +1,17 | |||||
|
1 | attenuation = numpy.array([[[-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
2 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
3 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
4 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
5 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
6 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
7 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
8 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25]], | |||
|
9 | [[21.25,21.25,21.25,21.25,21.25,21.25,21.25,21.25], | |||
|
10 | [15.25,15.25,15.25,15.25,15.25,15.25,15.25,15.25], | |||
|
11 | [09.25,09.25,09.25,09.25,09.25,09.25,09.25,09.25], | |||
|
12 | [03.25,03.25,03.25,03.25,03.25,03.25,03.25,03.25], | |||
|
13 | [-03.25,-03.25,-03.25,-03.25,-03.25,-03.25,-03.25,-03.25], | |||
|
14 | [-09.25,-09.25,-09.25,-09.25,-09.25,-09.25,-09.25,-09.25], | |||
|
15 | [-15.25,-15.25,-15.25,-15.25,-15.25,-15.25,-15.25,-15.25], | |||
|
16 | [-21.25,-21.25,-21.25,-21.25,-21.25,-21.25,-21.25,-21.25]]]) | |||
|
17 |
1 | NO CONTENT: new file 100644, binary diff hidden |
|
NO CONTENT: new file 100644, binary diff hidden |
@@ -0,0 +1,34 | |||||
|
1 | from numpy import * | |||
|
2 | from scipy import optimize | |||
|
3 | ||||
|
4 | def gaussian(height, center_x, center_y, width_x, width_y): | |||
|
5 | """Returns a gaussian function with the given parameters""" | |||
|
6 | width_x = float(width_x) | |||
|
7 | width_y = float(width_y) | |||
|
8 | return lambda x,y: height*exp( | |||
|
9 | -(((center_x-x)/width_x)**2+((center_y-y)/width_y)**2)/2) | |||
|
10 | ||||
|
11 | def moments(data): | |||
|
12 | """Returns (height, x, y, width_x, width_y) | |||
|
13 | the gaussian parameters of a 2D distribution by calculating its | |||
|
14 | moments """ | |||
|
15 | total = data.sum() | |||
|
16 | X, Y = indices(data.shape) | |||
|
17 | x = (X*data).sum()/total | |||
|
18 | y = (Y*data).sum()/total | |||
|
19 | col = data[:, int(y)] | |||
|
20 | width_x = sqrt(abs((arange(col.size)-y)**2*col).sum()/col.sum()) | |||
|
21 | row = data[int(x), :] | |||
|
22 | width_y = sqrt(abs((arange(row.size)-x)**2*row).sum()/row.sum()) | |||
|
23 | height = data.max() | |||
|
24 | return height, x, y, width_x, width_y | |||
|
25 | ||||
|
26 | def fitgaussian(data): | |||
|
27 | """Returns (height, x, y, width_x, width_y) | |||
|
28 | the gaussian parameters of a 2D distribution found by a fit""" | |||
|
29 | params = moments(data) | |||
|
30 | errorfunction = lambda p: ravel(gaussian(*p)(*indices(data.shape)) - | |||
|
31 | data) | |||
|
32 | p, success = optimize.leastsq(errorfunction, params) | |||
|
33 | return p | |||
|
34 |
1 | NO CONTENT: new file 100644, binary diff hidden |
|
NO CONTENT: new file 100644, binary diff hidden |
@@ -0,0 +1,196 | |||||
|
1 | g/h n m 1900.0 1905.0 1910.0 1915.0 1920.0 1925.0 1930.0 1935.0 1940.0 1945.0 1950.0 1955.0 1960.0 1965.0 1970.0 1975.0 1980.0 1985.0 1990.0 1995.0 2000.0 2005.0 SV | |||
|
2 | g 1 0 -31543 -31464 -31354 -31212 -31060 -30926 -30805 -30715 -30654 -30594 -30554 -30500 -30421 -30334 -30220 -30100 -29992 -29873 -29775 -29692 -29619.4 -29556.8 8.8 | |||
|
3 | g 1 1 -2298 -2298 -2297 -2306 -2317 -2318 -2316 -2306 -2292 -2285 -2250 -2215 -2169 -2119 -2068 -2013 -1956 -1905 -1848 -1784 -1728.2 -1671.8 10.8 | |||
|
4 | h 1 1 5922 5909 5898 5875 5845 5817 5808 5812 5821 5810 5815 5820 5791 5776 5737 5675 5604 5500 5406 5306 5186.1 5080.0 -21.3 | |||
|
5 | g 2 0 -677 -728 -769 -802 -839 -893 -951 -1018 -1106 -1244 -1341 -1440 -1555 -1662 -1781 -1902 -1997 -2072 -2131 -2200 -2267.7 -2340.5 -15.0 | |||
|
6 | g 2 1 2905 2928 2948 2956 2959 2969 2980 2984 2981 2990 2998 3003 3002 2997 3000 3010 3027 3044 3059 3070 3068.4 3047.0 -6.9 | |||
|
7 | h 2 1 -1061 -1086 -1128 -1191 -1259 -1334 -1424 -1520 -1614 -1702 -1810 -1898 -1967 -2016 -2047 -2067 -2129 -2197 -2279 -2366 -2481.6 -2594.9 -23.3 | |||
|
8 | g 2 2 924 1041 1176 1309 1407 1471 1517 1550 1566 1578 1576 1581 1590 1594 1611 1632 1663 1687 1686 1681 1670.9 1656.9 -1.0 | |||
|
9 | h 2 2 1121 1065 1000 917 823 728 644 586 528 477 381 291 206 114 25 -68 -200 -306 -373 -413 -458.0 -516.7 -14.0 | |||
|
10 | g 3 0 1022 1037 1058 1084 1111 1140 1172 1206 1240 1282 1297 1302 1302 1297 1287 1276 1281 1296 1314 1335 1339.6 1335.7 -0.3 | |||
|
11 | g 3 1 -1469 -1494 -1524 -1559 -1600 -1645 -1692 -1740 -1790 -1834 -1889 -1944 -1992 -2038 -2091 -2144 -2180 -2208 -2239 -2267 -2288.0 -2305.3 -3.1 | |||
|
12 | h 3 1 -330 -357 -389 -421 -445 -462 -480 -494 -499 -499 -476 -462 -414 -404 -366 -333 -336 -310 -284 -262 -227.6 -200.4 5.4 | |||
|
13 | g 3 2 1256 1239 1223 1212 1205 1202 1205 1215 1232 1255 1274 1288 1289 1292 1278 1260 1251 1247 1248 1249 1252.1 1246.8 -0.9 | |||
|
14 | h 3 2 3 34 62 84 103 119 133 146 163 186 206 216 224 240 251 262 271 284 293 302 293.4 269.3 -6.5 | |||
|
15 | g 3 3 572 635 705 778 839 881 907 918 916 913 896 882 878 856 838 830 833 829 802 759 714.5 674.4 -6.8 | |||
|
16 | h 3 3 523 480 425 360 293 229 166 101 43 -11 -46 -83 -130 -165 -196 -223 -252 -297 -352 -427 -491.1 -524.5 -2.0 | |||
|
17 | g 4 0 876 880 884 887 889 891 896 903 914 944 954 958 957 957 952 946 938 936 939 940 932.3 919.8 -2.5 | |||
|
18 | g 4 1 628 643 660 678 695 711 727 744 762 776 792 796 800 804 800 791 782 780 780 780 786.8 798.2 2.8 | |||
|
19 | h 4 1 195 203 211 218 220 216 205 188 169 144 136 133 135 148 167 191 212 232 247 262 272.6 281.4 2.0 | |||
|
20 | g 4 2 660 653 644 631 616 601 584 565 550 544 528 510 504 479 461 438 398 361 325 290 250.0 211.5 -7.1 | |||
|
21 | h 4 2 -69 -77 -90 -109 -134 -163 -195 -226 -252 -276 -278 -274 -278 -269 -266 -265 -257 -249 -240 -236 -231.9 -225.8 1.8 | |||
|
22 | g 4 3 -361 -380 -400 -416 -424 -426 -422 -415 -405 -421 -408 -397 -394 -390 -395 -405 -419 -424 -423 -418 -403.0 -379.5 5.9 | |||
|
23 | h 4 3 -210 -201 -189 -173 -153 -130 -109 -90 -72 -55 -37 -23 3 13 26 39 53 69 84 97 119.8 145.7 5.6 | |||
|
24 | g 4 4 134 146 160 178 199 217 234 249 265 304 303 290 269 252 234 216 199 170 141 122 111.3 100.2 -3.2 | |||
|
25 | h 4 4 -75 -65 -55 -51 -57 -70 -90 -114 -141 -178 -210 -230 -255 -269 -279 -288 -297 -297 -299 -306 -303.8 -304.7 0.0 | |||
|
26 | g 5 0 -184 -192 -201 -211 -221 -230 -237 -241 -241 -253 -240 -229 -222 -219 -216 -218 -218 -214 -214 -214 -218.8 -227.6 -2.6 | |||
|
27 | g 5 1 328 328 327 327 326 326 327 329 334 346 349 360 362 358 359 356 357 355 353 352 351.4 354.4 0.4 | |||
|
28 | h 5 1 -210 -193 -172 -148 -122 -96 -72 -51 -33 -12 3 15 16 19 26 31 46 47 46 46 43.8 42.7 0.1 | |||
|
29 | g 5 2 264 259 253 245 236 226 218 211 208 194 211 230 242 254 262 264 261 253 245 235 222.3 208.8 -3.0 | |||
|
30 | h 5 2 53 56 57 58 58 58 60 64 71 95 103 110 125 128 139 148 150 150 154 165 171.9 179.8 1.8 | |||
|
31 | g 5 3 5 -1 -9 -16 -23 -28 -32 -33 -33 -20 -20 -23 -26 -31 -42 -59 -74 -93 -109 -118 -130.4 -136.6 -1.2 | |||
|
32 | h 5 3 -33 -32 -33 -34 -38 -44 -53 -64 -75 -67 -87 -98 -117 -126 -139 -152 -151 -154 -153 -143 -133.1 -123.0 2.0 | |||
|
33 | g 5 4 -86 -93 -102 -111 -119 -125 -131 -136 -141 -142 -147 -152 -156 -157 -160 -159 -162 -164 -165 -166 -168.6 -168.3 0.2 | |||
|
34 | h 5 4 -124 -125 -126 -126 -125 -122 -118 -115 -113 -119 -122 -121 -114 -97 -91 -83 -78 -75 -69 -55 -39.3 -19.5 4.5 | |||
|
35 | g 5 5 -16 -26 -38 -51 -62 -69 -74 -76 -76 -82 -76 -69 -63 -62 -56 -49 -48 -46 -36 -17 -12.9 -14.1 -0.6 | |||
|
36 | h 5 5 3 11 21 32 43 51 58 64 69 82 80 78 81 81 83 88 92 95 97 107 106.3 103.6 -1.0 | |||
|
37 | g 6 0 63 62 62 61 61 61 60 59 57 59 54 47 46 45 43 45 48 53 61 68 72.3 72.9 -0.8 | |||
|
38 | g 6 1 61 60 58 57 55 54 53 53 54 57 57 57 58 61 64 66 66 65 65 67 68.2 69.6 0.2 | |||
|
39 | h 6 1 -9 -7 -5 -2 0 3 4 4 4 6 -1 -9 -10 -11 -12 -13 -15 -16 -16 -17 -17.4 -20.2 -0.4 | |||
|
40 | g 6 2 -11 -11 -11 -10 -10 -9 -9 -8 -7 6 4 3 1 8 15 28 42 51 59 68 74.2 76.6 -0.2 | |||
|
41 | h 6 2 83 86 89 93 96 99 102 104 105 100 99 96 99 100 100 99 93 88 82 72 63.7 54.7 -1.9 | |||
|
42 | g 6 3 -217 -221 -224 -228 -233 -238 -242 -246 -249 -246 -247 -247 -237 -228 -212 -198 -192 -185 -178 -170 -160.9 -151.1 2.1 | |||
|
43 | h 6 3 2 4 5 8 11 14 19 25 33 16 33 48 60 68 72 75 71 69 69 67 65.1 63.7 -0.4 | |||
|
44 | g 6 4 -58 -57 -54 -51 -46 -40 -32 -25 -18 -25 -16 -8 -1 4 2 1 4 4 3 -1 -5.9 -15.0 -2.1 | |||
|
45 | h 6 4 -35 -32 -29 -26 -22 -18 -16 -15 -15 -9 -12 -16 -20 -32 -37 -41 -43 -48 -52 -58 -61.2 -63.4 -0.4 | |||
|
46 | g 6 5 59 57 54 49 44 39 32 25 18 21 12 7 -2 1 3 6 14 16 18 19 16.9 14.7 -0.4 | |||
|
47 | h 6 5 36 32 28 23 18 13 8 4 0 -16 -12 -12 -11 -8 -6 -4 -2 -1 1 1 0.7 0.0 -0.2 | |||
|
48 | g 6 6 -90 -92 -95 -98 -101 -103 -104 -106 -107 -104 -105 -107 -113 -111 -112 -111 -108 -102 -96 -93 -90.4 -86.4 1.3 | |||
|
49 | h 6 6 -69 -67 -65 -62 -57 -52 -46 -40 -33 -39 -30 -24 -17 -7 1 11 17 21 24 36 43.8 50.3 0.9 | |||
|
50 | g 7 0 70 70 71 72 73 73 74 74 74 70 65 65 67 75 72 71 72 74 77 77 79.0 79.8 -0.4 | |||
|
51 | g 7 1 -55 -54 -54 -54 -54 -54 -54 -53 -53 -40 -55 -56 -56 -57 -57 -56 -59 -62 -64 -72 -74.0 -74.4 0.0 | |||
|
52 | h 7 1 -45 -46 -47 -48 -49 -50 -51 -52 -52 -45 -35 -50 -55 -61 -70 -77 -82 -83 -80 -69 -64.6 -61.4 0.8 | |||
|
53 | g 7 2 0 0 1 2 2 3 4 4 4 0 2 2 5 4 1 1 2 3 2 1 0.0 -1.4 -0.2 | |||
|
54 | h 7 2 -13 -14 -14 -14 -14 -14 -15 -17 -18 -18 -17 -24 -28 -27 -27 -26 -27 -27 -26 -25 -24.2 -22.5 0.4 | |||
|
55 | g 7 3 34 33 32 31 29 27 25 23 20 0 1 10 15 13 14 16 21 24 26 28 33.3 38.6 1.1 | |||
|
56 | h 7 3 -10 -11 -12 -12 -13 -14 -14 -14 -14 2 0 -4 -6 -2 -4 -5 -5 -2 0 4 6.2 6.9 0.1 | |||
|
57 | g 7 4 -41 -41 -40 -38 -37 -35 -34 -33 -31 -29 -40 -32 -32 -26 -22 -14 -12 -6 -1 5 9.1 12.3 0.6 | |||
|
58 | h 7 4 -1 0 1 2 4 5 6 7 7 6 10 8 7 6 8 10 16 20 21 24 24.0 25.4 0.2 | |||
|
59 | g 7 5 -21 -20 -19 -18 -16 -14 -12 -11 -9 -10 -7 -11 -7 -6 -2 0 1 4 5 4 6.9 9.4 0.4 | |||
|
60 | h 7 5 28 28 28 28 28 29 29 29 29 28 36 28 23 26 23 22 18 17 17 17 14.8 10.9 -0.9 | |||
|
61 | g 7 6 18 18 18 19 19 19 18 18 17 15 5 9 17 13 13 12 11 10 9 8 7.3 5.5 -0.5 | |||
|
62 | h 7 6 -12 -12 -13 -15 -16 -17 -18 -19 -20 -17 -18 -20 -18 -23 -23 -23 -23 -23 -23 -24 -25.4 -26.4 -0.3 | |||
|
63 | g 7 7 6 6 6 6 6 6 6 6 5 29 19 18 8 1 -2 -5 -2 0 0 -2 -1.2 2.0 0.9 | |||
|
64 | h 7 7 -22 -22 -22 -22 -22 -21 -20 -19 -19 -22 -16 -18 -17 -12 -11 -12 -10 -7 -4 -6 -5.8 -4.8 0.3 | |||
|
65 | g 8 0 11 11 11 11 11 11 11 11 11 13 22 11 15 13 14 14 18 21 23 25 24.4 24.8 -0.2 | |||
|
66 | g 8 1 8 8 8 8 7 7 7 7 7 7 15 9 6 5 6 6 6 6 5 6 6.6 7.7 0.2 | |||
|
67 | h 8 1 8 8 8 8 8 8 8 8 8 12 5 10 11 7 7 6 7 8 10 11 11.9 11.2 -0.2 | |||
|
68 | g 8 2 -4 -4 -4 -4 -3 -3 -3 -3 -3 -8 -4 -6 -4 -4 -2 -1 0 0 -1 -6 -9.2 -11.4 -0.2 | |||
|
69 | h 8 2 -14 -15 -15 -15 -15 -15 -15 -15 -14 -21 -22 -15 -14 -12 -15 -16 -18 -19 -19 -21 -21.5 -21.0 0.2 | |||
|
70 | g 8 3 -9 -9 -9 -9 -9 -9 -9 -9 -10 -5 -1 -14 -11 -14 -13 -12 -11 -11 -10 -9 -7.9 -6.8 0.2 | |||
|
71 | h 8 3 7 7 6 6 6 6 5 5 5 -12 0 5 7 9 6 4 4 5 6 8 8.5 9.7 0.2 | |||
|
72 | g 8 4 1 1 1 2 2 2 2 1 1 9 11 6 2 0 -3 -8 -7 -9 -12 -14 -16.6 -18.0 -0.2 | |||
|
73 | h 8 4 -13 -13 -13 -13 -14 -14 -14 -15 -15 -7 -21 -23 -18 -16 -17 -19 -22 -23 -22 -23 -21.5 -19.8 0.4 | |||
|
74 | g 8 5 2 2 2 3 4 4 5 6 6 7 15 10 10 8 5 4 4 4 3 9 9.1 10.0 0.2 | |||
|
75 | h 8 5 5 5 5 5 5 5 5 5 5 2 -8 3 4 4 6 6 9 11 12 15 15.5 16.1 0.2 | |||
|
76 | g 8 6 -9 -8 -8 -8 -7 -7 -6 -6 -5 -10 -13 -7 -5 -1 0 0 3 4 4 6 7.0 9.4 0.5 | |||
|
77 | h 8 6 16 16 16 16 17 17 18 18 19 18 17 23 23 24 21 18 16 14 12 11 8.9 7.7 -0.3 | |||
|
78 | g 8 7 5 5 5 6 6 7 8 8 9 7 5 6 10 11 11 10 6 4 2 -5 -7.9 -11.4 -0.7 | |||
|
79 | h 8 7 -5 -5 -5 -5 -5 -5 -5 -5 -5 3 -4 -4 1 -3 -6 -10 -13 -15 -16 -16 -14.9 -12.8 0.5 | |||
|
80 | g 8 8 8 8 8 8 8 8 8 7 7 2 -1 9 8 4 3 1 -1 -4 -6 -7 -7.0 -5.0 0.5 | |||
|
81 | h 8 8 -18 -18 -18 -18 -19 -19 -19 -19 -19 -11 -17 -13 -20 -17 -16 -17 -15 -11 -10 -4 -2.1 -0.1 0.4 | |||
|
82 | g 9 0 8 8 8 8 8 8 8 8 8 5 3 4 4 8 8 7 5 5 4 4 5.0 5.6 | |||
|
83 | g 9 1 10 10 10 10 10 10 10 10 10 -21 -7 9 6 10 10 10 10 10 9 9 9.4 9.8 | |||
|
84 | h 9 1 -20 -20 -20 -20 -20 -20 -20 -20 -21 -27 -24 -11 -18 -22 -21 -21 -21 -21 -20 -20 -19.7 -20.1 | |||
|
85 | g 9 2 1 1 1 1 1 1 1 1 1 1 -1 -4 0 2 2 2 1 1 1 3 3.0 3.6 | |||
|
86 | h 9 2 14 14 14 14 14 14 14 15 15 17 19 12 12 15 16 16 16 15 15 15 13.4 12.9 | |||
|
87 | g 9 3 -11 -11 -11 -11 -11 -11 -12 -12 -12 -11 -25 -5 -9 -13 -12 -12 -12 -12 -12 -10 -8.4 -7.0 | |||
|
88 | h 9 3 5 5 5 5 5 5 5 5 5 29 12 7 2 7 6 7 9 9 11 12 12.5 12.7 | |||
|
89 | g 9 4 12 12 12 12 12 12 12 11 11 3 10 2 1 10 10 10 9 9 9 8 6.3 5.0 | |||
|
90 | h 9 4 -3 -3 -3 -3 -3 -3 -3 -3 -3 -9 2 6 0 -4 -4 -4 -5 -6 -7 -6 -6.2 -6.7 | |||
|
91 | g 9 5 1 1 1 1 1 1 1 1 1 16 5 4 4 -1 -1 -1 -3 -3 -4 -8 -8.9 -10.8 | |||
|
92 | h 9 5 -2 -2 -2 -2 -2 -2 -2 -3 -3 4 2 -2 -3 -5 -5 -5 -6 -6 -7 -8 -8.4 -8.1 | |||
|
93 | g 9 6 -2 -2 -2 -2 -2 -2 -2 -2 -2 -3 -5 1 -1 -1 0 -1 -1 -1 -2 -1 -1.5 -1.3 | |||
|
94 | h 9 6 8 8 8 8 9 9 9 9 9 9 8 10 9 10 10 10 9 9 9 8 8.4 8.1 | |||
|
95 | g 9 7 2 2 2 2 2 2 3 3 3 -4 -2 2 -2 5 3 4 7 7 7 10 9.3 8.7 | |||
|
96 | h 9 7 10 10 10 10 10 10 10 11 11 6 8 7 8 10 11 11 10 9 8 5 3.8 2.9 | |||
|
97 | g 9 8 -1 0 0 0 0 0 0 0 1 -3 3 2 3 1 1 1 2 1 1 -2 -4.3 -6.7 | |||
|
98 | h 9 8 -2 -2 -2 -2 -2 -2 -2 -2 -2 1 -11 -6 0 -4 -2 -3 -6 -7 -7 -8 -8.2 -7.9 | |||
|
99 | g 9 9 -1 -1 -1 -1 -1 -1 -2 -2 -2 -4 8 5 -1 -2 -1 -2 -5 -5 -6 -8 -8.2 -9.2 | |||
|
100 | h 9 9 2 2 2 2 2 2 2 2 2 8 -7 5 5 1 1 1 2 2 2 3 4.8 5.9 | |||
|
101 | g 10 0 -3 -3 -3 -3 -3 -3 -3 -3 -3 -3 -8 -3 1 -2 -3 -3 -4 -4 -3 -3 -2.6 -2.2 | |||
|
102 | g 10 1 -4 -4 -4 -4 -4 -4 -4 -4 -4 11 4 -5 -3 -3 -3 -3 -4 -4 -4 -6 -6.0 -6.3 | |||
|
103 | h 10 1 2 2 2 2 2 2 2 2 2 5 13 -4 4 2 1 1 1 1 2 1 1.7 2.4 | |||
|
104 | g 10 2 2 2 2 2 2 2 2 2 2 1 -1 -1 4 2 2 2 2 3 2 2 1.7 1.6 | |||
|
105 | h 10 2 1 1 1 1 1 1 1 1 1 1 -2 0 1 1 1 1 0 0 1 0 0.0 0.2 | |||
|
106 | g 10 3 -5 -5 -5 -5 -5 -5 -5 -5 -5 2 13 2 0 -5 -5 -5 -5 -5 -5 -4 -3.1 -2.5 | |||
|
107 | h 10 3 2 2 2 2 2 2 2 2 2 -20 -10 -8 0 2 3 3 3 3 3 4 4.0 4.4 | |||
|
108 | g 10 4 -2 -2 -2 -2 -2 -2 -2 -2 -2 -5 -4 -3 -1 -2 -1 -2 -2 -2 -2 -1 -0.5 -0.1 | |||
|
109 | h 10 4 6 6 6 6 6 6 6 6 6 -1 2 -2 2 6 4 4 6 6 6 5 4.9 4.7 | |||
|
110 | g 10 5 6 6 6 6 6 6 6 6 6 -1 4 7 4 4 6 5 5 5 4 4 3.7 3.0 | |||
|
111 | h 10 5 -4 -4 -4 -4 -4 -4 -4 -4 -4 -6 -3 -4 -5 -4 -4 -4 -4 -4 -4 -5 -5.9 -6.5 | |||
|
112 | g 10 6 4 4 4 4 4 4 4 4 4 8 12 4 6 4 4 4 3 3 3 2 1.0 0.3 | |||
|
113 | h 10 6 0 0 0 0 0 0 0 0 0 6 6 1 1 0 0 -1 0 0 0 -1 -1.2 -1.0 | |||
|
114 | g 10 7 0 0 0 0 0 0 0 0 0 -1 3 -2 1 0 1 1 1 1 1 2 2.0 2.1 | |||
|
115 | h 10 7 -2 -2 -2 -2 -2 -2 -2 -1 -1 -4 -3 -3 -1 -2 -1 -1 -1 -1 -2 -2 -2.9 -3.4 | |||
|
116 | g 10 8 2 2 2 1 1 1 1 2 2 -3 2 6 -1 2 0 0 2 2 3 5 4.2 3.9 | |||
|
117 | h 10 8 4 4 4 4 4 4 4 4 4 -2 6 7 6 3 3 3 4 4 3 1 0.2 -0.9 | |||
|
118 | g 10 9 2 2 2 2 3 3 3 3 3 5 10 -2 2 2 3 3 3 3 3 1 0.3 -0.1 | |||
|
119 | h 10 9 0 0 0 0 0 0 0 0 0 0 11 -1 0 0 1 1 0 0 -1 -2 -2.2 -2.3 | |||
|
120 | g 10 10 0 0 0 0 0 0 0 0 0 -2 3 0 0 0 -1 -1 0 0 0 0 -1.1 -2.2 | |||
|
121 | h 10 10 -6 -6 -6 -6 -6 -6 -6 -6 -6 -2 8 -3 -7 -6 -4 -5 -6 -6 -6 -7 -7.4 -8.0 | |||
|
122 | g 11 0 2.7 2.9 | |||
|
123 | g 11 1 -1.7 -1.6 | |||
|
124 | h 11 1 0.1 0.3 | |||
|
125 | g 11 2 -1.9 -1.7 | |||
|
126 | h 11 2 1.3 1.4 | |||
|
127 | g 11 3 1.5 1.5 | |||
|
128 | h 11 3 -0.9 -0.7 | |||
|
129 | g 11 4 -0.1 -0.2 | |||
|
130 | h 11 4 -2.6 -2.4 | |||
|
131 | g 11 5 0.1 0.2 | |||
|
132 | h 11 5 0.9 0.9 | |||
|
133 | g 11 6 -0.7 -0.7 | |||
|
134 | h 11 6 -0.7 -0.6 | |||
|
135 | g 11 7 0.7 0.5 | |||
|
136 | h 11 7 -2.8 -2.7 | |||
|
137 | g 11 8 1.7 1.8 | |||
|
138 | h 11 8 -0.9 -1.0 | |||
|
139 | g 11 9 0.1 0.1 | |||
|
140 | h 11 9 -1.2 -1.5 | |||
|
141 | g 11 10 1.2 1.0 | |||
|
142 | h 11 10 -1.9 -2.0 | |||
|
143 | g 11 11 4.0 4.1 | |||
|
144 | h 11 11 -0.9 -1.4 | |||
|
145 | g 12 0 -2.2 -2.2 | |||
|
146 | g 12 1 -0.3 -0.3 | |||
|
147 | h 12 1 -0.4 -0.5 | |||
|
148 | g 12 2 0.2 0.3 | |||
|
149 | h 12 2 0.3 0.3 | |||
|
150 | g 12 3 0.9 0.9 | |||
|
151 | h 12 3 2.5 2.3 | |||
|
152 | g 12 4 -0.2 -0.4 | |||
|
153 | h 12 4 -2.6 -2.7 | |||
|
154 | g 12 5 0.9 1.0 | |||
|
155 | h 12 5 0.7 0.6 | |||
|
156 | g 12 6 -0.5 -0.4 | |||
|
157 | h 12 6 0.3 0.4 | |||
|
158 | g 12 7 0.3 0.5 | |||
|
159 | h 12 7 0.0 0.0 | |||
|
160 | g 12 8 -0.3 -0.3 | |||
|
161 | h 12 8 0.0 0.0 | |||
|
162 | g 12 9 -0.4 -0.4 | |||
|
163 | h 12 9 0.3 0.3 | |||
|
164 | g 12 10 -0.1 0.0 | |||
|
165 | h 12 10 -0.9 -0.8 | |||
|
166 | g 12 11 -0.2 -0.4 | |||
|
167 | h 12 11 -0.4 -0.4 | |||
|
168 | g 12 12 -0.4 0.0 | |||
|
169 | h 12 12 0.8 1.0 | |||
|
170 | g 13 0 -0.2 -0.2 | |||
|
171 | g 13 1 -0.9 -0.9 | |||
|
172 | h 13 1 -0.9 -0.7 | |||
|
173 | g 13 2 0.3 0.3 | |||
|
174 | h 13 2 0.2 0.3 | |||
|
175 | g 13 3 0.1 0.3 | |||
|
176 | h 13 3 1.8 1.7 | |||
|
177 | g 13 4 -0.4 -0.4 | |||
|
178 | h 13 4 -0.4 -0.5 | |||
|
179 | g 13 5 1.3 1.2 | |||
|
180 | h 13 5 -1.0 -1.0 | |||
|
181 | g 13 6 -0.4 -0.4 | |||
|
182 | h 13 6 -0.1 0.0 | |||
|
183 | g 13 7 0.7 0.7 | |||
|
184 | h 13 7 0.7 0.7 | |||
|
185 | g 13 8 -0.4 -0.3 | |||
|
186 | h 13 8 0.3 0.2 | |||
|
187 | g 13 9 0.3 0.4 | |||
|
188 | h 13 9 0.6 0.6 | |||
|
189 | g 13 10 -0.1 -0.1 | |||
|
190 | h 13 10 0.3 0.4 | |||
|
191 | g 13 11 0.4 0.4 | |||
|
192 | h 13 11 -0.2 -0.2 | |||
|
193 | g 13 12 0.0 -0.1 | |||
|
194 | h 13 12 -0.5 -0.5 | |||
|
195 | g 13 13 0.1 -0.3 | |||
|
196 | h 13 13 -0.9 -1.0 |
@@ -0,0 +1,196 | |||||
|
1 | g/h n m 1900.0 1905.0 1910.0 1915.0 1920.0 1925.0 1930.0 1935.0 1940.0 1945.0 1950.0 1955.0 1960.0 1965.0 1970.0 1975.0 1980.0 1985.0 1990.0 1995.0 2000.0 2005.0 2010.0 SV | |||
|
2 | g 1 0 -31543 -31464 -31354 -31212 -31060 -30926 -30805 -30715 -30654 -30594 -30554 -30500 -30421 -30334 -30220 -30100 -29992 -29873 -29775 -29692 -29619.4 -29554.63 -29496.5 11.4 | |||
|
3 | g 1 1 -2298 -2298 -2297 -2306 -2317 -2318 -2316 -2306 -2292 -2285 -2250 -2215 -2169 -2119 -2068 -2013 -1956 -1905 -1848 -1784 -1728.2 -1669.05 -1585.9 16.7 | |||
|
4 | h 1 1 5922 5909 5898 5875 5845 5817 5808 5812 5821 5810 5815 5820 5791 5776 5737 5675 5604 5500 5406 5306 5186.1 5077.99 4945.1 -28.8 | |||
|
5 | g 2 0 -677 -728 -769 -802 -839 -893 -951 -1018 -1106 -1244 -1341 -1440 -1555 -1662 -1781 -1902 -1997 -2072 -2131 -2200 -2267.7 -2337.24 -2396.6 -11.3 | |||
|
6 | g 2 1 2905 2928 2948 2956 2959 2969 2980 2984 2981 2990 2998 3003 3002 2997 3000 3010 3027 3044 3059 3070 3068.4 3047.69 3026.0 -3.9 | |||
|
7 | h 2 1 -1061 -1086 -1128 -1191 -1259 -1334 -1424 -1520 -1614 -1702 -1810 -1898 -1967 -2016 -2047 -2067 -2129 -2197 -2279 -2366 -2481.6 -2594.50 -2707.7 -23.0 | |||
|
8 | g 2 2 924 1041 1176 1309 1407 1471 1517 1550 1566 1578 1576 1581 1590 1594 1611 1632 1663 1687 1686 1681 1670.9 1657.76 1668.6 2.7 | |||
|
9 | h 2 2 1121 1065 1000 917 823 728 644 586 528 477 381 291 206 114 25 -68 -200 -306 -373 -413 -458.0 -515.43 -575.4 -12.9 | |||
|
10 | g 3 0 1022 1037 1058 1084 1111 1140 1172 1206 1240 1282 1297 1302 1302 1297 1287 1276 1281 1296 1314 1335 1339.6 1336.30 1339.7 1.3 | |||
|
11 | g 3 1 -1469 -1494 -1524 -1559 -1600 -1645 -1692 -1740 -1790 -1834 -1889 -1944 -1992 -2038 -2091 -2144 -2180 -2208 -2239 -2267 -2288.0 -2305.83 -2326.3 -3.9 | |||
|
12 | h 3 1 -330 -357 -389 -421 -445 -462 -480 -494 -499 -499 -476 -462 -414 -404 -366 -333 -336 -310 -284 -262 -227.6 -198.86 -160.5 8.6 | |||
|
13 | g 3 2 1256 1239 1223 1212 1205 1202 1205 1215 1232 1255 1274 1288 1289 1292 1278 1260 1251 1247 1248 1249 1252.1 1246.39 1231.7 -2.9 | |||
|
14 | h 3 2 3 34 62 84 103 119 133 146 163 186 206 216 224 240 251 262 271 284 293 302 293.4 269.72 251.7 -2.9 | |||
|
15 | g 3 3 572 635 705 778 839 881 907 918 916 913 896 882 878 856 838 830 833 829 802 759 714.5 672.51 634.2 -8.1 | |||
|
16 | h 3 3 523 480 425 360 293 229 166 101 43 -11 -46 -83 -130 -165 -196 -223 -252 -297 -352 -427 -491.1 -524.72 -536.8 -2.1 | |||
|
17 | g 4 0 876 880 884 887 889 891 896 903 914 944 954 958 957 957 952 946 938 936 939 940 932.3 920.55 912.6 -1.4 | |||
|
18 | g 4 1 628 643 660 678 695 711 727 744 762 776 792 796 800 804 800 791 782 780 780 780 786.8 797.96 809.0 2.0 | |||
|
19 | h 4 1 195 203 211 218 220 216 205 188 169 144 136 133 135 148 167 191 212 232 247 262 272.6 282.07 286.4 0.4 | |||
|
20 | g 4 2 660 653 644 631 616 601 584 565 550 544 528 510 504 479 461 438 398 361 325 290 250.0 210.65 166.6 -8.9 | |||
|
21 | h 4 2 -69 -77 -90 -109 -134 -163 -195 -226 -252 -276 -278 -274 -278 -269 -266 -265 -257 -249 -240 -236 -231.9 -225.23 -211.2 3.2 | |||
|
22 | g 4 3 -361 -380 -400 -416 -424 -426 -422 -415 -405 -421 -408 -397 -394 -390 -395 -405 -419 -424 -423 -418 -403.0 -379.86 -357.1 4.4 | |||
|
23 | h 4 3 -210 -201 -189 -173 -153 -130 -109 -90 -72 -55 -37 -23 3 13 26 39 53 69 84 97 119.8 145.15 164.4 3.6 | |||
|
24 | g 4 4 134 146 160 178 199 217 234 249 265 304 303 290 269 252 234 216 199 170 141 122 111.3 100.00 89.7 -2.3 | |||
|
25 | h 4 4 -75 -65 -55 -51 -57 -70 -90 -114 -141 -178 -210 -230 -255 -269 -279 -288 -297 -297 -299 -306 -303.8 -305.36 -309.2 -0.8 | |||
|
26 | g 5 0 -184 -192 -201 -211 -221 -230 -237 -241 -241 -253 -240 -229 -222 -219 -216 -218 -218 -214 -214 -214 -218.8 -227.00 -231.1 -0.5 | |||
|
27 | g 5 1 328 328 327 327 326 326 327 329 334 346 349 360 362 358 359 356 357 355 353 352 351.4 354.41 357.2 0.5 | |||
|
28 | h 5 1 -210 -193 -172 -148 -122 -96 -72 -51 -33 -12 3 15 16 19 26 31 46 47 46 46 43.8 42.72 44.7 0.5 | |||
|
29 | g 5 2 264 259 253 245 236 226 218 211 208 194 211 230 242 254 262 264 261 253 245 235 222.3 208.95 200.3 -1.5 | |||
|
30 | h 5 2 53 56 57 58 58 58 60 64 71 95 103 110 125 128 139 148 150 150 154 165 171.9 180.25 188.9 1.5 | |||
|
31 | g 5 3 5 -1 -9 -16 -23 -28 -32 -33 -33 -20 -20 -23 -26 -31 -42 -59 -74 -93 -109 -118 -130.4 -136.54 -141.2 -0.7 | |||
|
32 | h 5 3 -33 -32 -33 -34 -38 -44 -53 -64 -75 -67 -87 -98 -117 -126 -139 -152 -151 -154 -153 -143 -133.1 -123.45 -118.1 0.9 | |||
|
33 | g 5 4 -86 -93 -102 -111 -119 -125 -131 -136 -141 -142 -147 -152 -156 -157 -160 -159 -162 -164 -165 -166 -168.6 -168.05 -163.1 1.3 | |||
|
34 | h 5 4 -124 -125 -126 -126 -125 -122 -118 -115 -113 -119 -122 -121 -114 -97 -91 -83 -78 -75 -69 -55 -39.3 -19.57 0.1 3.7 | |||
|
35 | g 5 5 -16 -26 -38 -51 -62 -69 -74 -76 -76 -82 -76 -69 -63 -62 -56 -49 -48 -46 -36 -17 -12.9 -13.55 -7.7 1.4 | |||
|
36 | h 5 5 3 11 21 32 43 51 58 64 69 82 80 78 81 81 83 88 92 95 97 107 106.3 103.85 100.9 -0.6 | |||
|
37 | g 6 0 63 62 62 61 61 61 60 59 57 59 54 47 46 45 43 45 48 53 61 68 72.3 73.60 72.8 -0.3 | |||
|
38 | g 6 1 61 60 58 57 55 54 53 53 54 57 57 57 58 61 64 66 66 65 65 67 68.2 69.56 68.6 -0.3 | |||
|
39 | h 6 1 -9 -7 -5 -2 0 3 4 4 4 6 -1 -9 -10 -11 -12 -13 -15 -16 -16 -17 -17.4 -20.33 -20.8 -0.1 | |||
|
40 | g 6 2 -11 -11 -11 -10 -10 -9 -9 -8 -7 6 4 3 1 8 15 28 42 51 59 68 74.2 76.74 76.0 -0.3 | |||
|
41 | h 6 2 83 86 89 93 96 99 102 104 105 100 99 96 99 100 100 99 93 88 82 72 63.7 54.75 44.2 -2.1 | |||
|
42 | g 6 3 -217 -221 -224 -228 -233 -238 -242 -246 -249 -246 -247 -247 -237 -228 -212 -198 -192 -185 -178 -170 -160.9 -151.34 -141.4 1.9 | |||
|
43 | h 6 3 2 4 5 8 11 14 19 25 33 16 33 48 60 68 72 75 71 69 69 67 65.1 63.63 61.5 -0.4 | |||
|
44 | g 6 4 -58 -57 -54 -51 -46 -40 -32 -25 -18 -25 -16 -8 -1 4 2 1 4 4 3 -1 -5.9 -14.58 -22.9 -1.6 | |||
|
45 | h 6 4 -35 -32 -29 -26 -22 -18 -16 -15 -15 -9 -12 -16 -20 -32 -37 -41 -43 -48 -52 -58 -61.2 -63.53 -66.3 -0.5 | |||
|
46 | g 6 5 59 57 54 49 44 39 32 25 18 21 12 7 -2 1 3 6 14 16 18 19 16.9 14.58 13.1 -0.2 | |||
|
47 | h 6 5 36 32 28 23 18 13 8 4 0 -16 -12 -12 -11 -8 -6 -4 -2 -1 1 1 0.7 0.24 3.1 0.8 | |||
|
48 | g 6 6 -90 -92 -95 -98 -101 -103 -104 -106 -107 -104 -105 -107 -113 -111 -112 -111 -108 -102 -96 -93 -90.4 -86.36 -77.9 1.8 | |||
|
49 | h 6 6 -69 -67 -65 -62 -57 -52 -46 -40 -33 -39 -30 -24 -17 -7 1 11 17 21 24 36 43.8 50.94 54.9 0.5 | |||
|
50 | g 7 0 70 70 71 72 73 73 74 74 74 70 65 65 67 75 72 71 72 74 77 77 79.0 79.88 80.4 0.2 | |||
|
51 | g 7 1 -55 -54 -54 -54 -54 -54 -54 -53 -53 -40 -55 -56 -56 -57 -57 -56 -59 -62 -64 -72 -74.0 -74.46 -75.0 -0.1 | |||
|
52 | h 7 1 -45 -46 -47 -48 -49 -50 -51 -52 -52 -45 -35 -50 -55 -61 -70 -77 -82 -83 -80 -69 -64.6 -61.14 -57.8 0.6 | |||
|
53 | g 7 2 0 0 1 2 2 3 4 4 4 0 2 2 5 4 1 1 2 3 2 1 0.0 -1.65 -4.7 -0.6 | |||
|
54 | h 7 2 -13 -14 -14 -14 -14 -14 -15 -17 -18 -18 -17 -24 -28 -27 -27 -26 -27 -27 -26 -25 -24.2 -22.57 -21.2 0.3 | |||
|
55 | g 7 3 34 33 32 31 29 27 25 23 20 0 1 10 15 13 14 16 21 24 26 28 33.3 38.73 45.3 1.4 | |||
|
56 | h 7 3 -10 -11 -12 -12 -13 -14 -14 -14 -14 2 0 -4 -6 -2 -4 -5 -5 -2 0 4 6.2 6.82 6.6 -0.2 | |||
|
57 | g 7 4 -41 -41 -40 -38 -37 -35 -34 -33 -31 -29 -40 -32 -32 -26 -22 -14 -12 -6 -1 5 9.1 12.30 14.0 0.3 | |||
|
58 | h 7 4 -1 0 1 2 4 5 6 7 7 6 10 8 7 6 8 10 16 20 21 24 24.0 25.35 24.9 -0.1 | |||
|
59 | g 7 5 -21 -20 -19 -18 -16 -14 -12 -11 -9 -10 -7 -11 -7 -6 -2 0 1 4 5 4 6.9 9.37 10.4 0.1 | |||
|
60 | h 7 5 28 28 28 28 28 29 29 29 29 28 36 28 23 26 23 22 18 17 17 17 14.8 10.93 7.0 -0.8 | |||
|
61 | g 7 6 18 18 18 19 19 19 18 18 17 15 5 9 17 13 13 12 11 10 9 8 7.3 5.42 1.6 -0.8 | |||
|
62 | h 7 6 -12 -12 -13 -15 -16 -17 -18 -19 -20 -17 -18 -20 -18 -23 -23 -23 -23 -23 -23 -24 -25.4 -26.32 -27.7 -0.3 | |||
|
63 | g 7 7 6 6 6 6 6 6 6 6 5 29 19 18 8 1 -2 -5 -2 0 0 -2 -1.2 1.94 4.9 0.4 | |||
|
64 | h 7 7 -22 -22 -22 -22 -22 -21 -20 -19 -19 -22 -16 -18 -17 -12 -11 -12 -10 -7 -4 -6 -5.8 -4.64 -3.4 0.2 | |||
|
65 | g 8 0 11 11 11 11 11 11 11 11 11 13 22 11 15 13 14 14 18 21 23 25 24.4 24.80 24.3 -0.1 | |||
|
66 | g 8 1 8 8 8 8 7 7 7 7 7 7 15 9 6 5 6 6 6 6 5 6 6.6 7.62 8.2 0.1 | |||
|
67 | h 8 1 8 8 8 8 8 8 8 8 8 12 5 10 11 7 7 6 7 8 10 11 11.9 11.20 10.9 0.0 | |||
|
68 | g 8 2 -4 -4 -4 -4 -3 -3 -3 -3 -3 -8 -4 -6 -4 -4 -2 -1 0 0 -1 -6 -9.2 -11.73 -14.5 -0.5 | |||
|
69 | h 8 2 -14 -15 -15 -15 -15 -15 -15 -15 -14 -21 -22 -15 -14 -12 -15 -16 -18 -19 -19 -21 -21.5 -20.88 -20.0 0.2 | |||
|
70 | g 8 3 -9 -9 -9 -9 -9 -9 -9 -9 -10 -5 -1 -14 -11 -14 -13 -12 -11 -11 -10 -9 -7.9 -6.88 -5.7 0.3 | |||
|
71 | h 8 3 7 7 6 6 6 6 5 5 5 -12 0 5 7 9 6 4 4 5 6 8 8.5 9.83 11.9 0.5 | |||
|
72 | g 8 4 1 1 1 2 2 2 2 1 1 9 11 6 2 0 -3 -8 -7 -9 -12 -14 -16.6 -18.11 -19.3 -0.3 | |||
|
73 | h 8 4 -13 -13 -13 -13 -14 -14 -14 -15 -15 -7 -21 -23 -18 -16 -17 -19 -22 -23 -22 -23 -21.5 -19.71 -17.4 0.4 | |||
|
74 | g 8 5 2 2 2 3 4 4 5 6 6 7 15 10 10 8 5 4 4 4 3 9 9.1 10.17 11.6 0.3 | |||
|
75 | h 8 5 5 5 5 5 5 5 5 5 5 2 -8 3 4 4 6 6 9 11 12 15 15.5 16.22 16.7 0.1 | |||
|
76 | g 8 6 -9 -8 -8 -8 -7 -7 -6 -6 -5 -10 -13 -7 -5 -1 0 0 3 4 4 6 7.0 9.36 10.9 0.2 | |||
|
77 | h 8 6 16 16 16 16 17 17 18 18 19 18 17 23 23 24 21 18 16 14 12 11 8.9 7.61 7.1 -0.1 | |||
|
78 | g 8 7 5 5 5 6 6 7 8 8 9 7 5 6 10 11 11 10 6 4 2 -5 -7.9 -11.25 -14.1 -0.5 | |||
|
79 | h 8 7 -5 -5 -5 -5 -5 -5 -5 -5 -5 3 -4 -4 1 -3 -6 -10 -13 -15 -16 -16 -14.9 -12.76 -10.8 0.4 | |||
|
80 | g 8 8 8 8 8 8 8 8 8 7 7 2 -1 9 8 4 3 1 -1 -4 -6 -7 -7.0 -4.87 -3.7 0.2 | |||
|
81 | h 8 8 -18 -18 -18 -18 -19 -19 -19 -19 -19 -11 -17 -13 -20 -17 -16 -17 -15 -11 -10 -4 -2.1 -0.06 1.7 0.4 | |||
|
82 | g 9 0 8 8 8 8 8 8 8 8 8 5 3 4 4 8 8 7 5 5 4 4 5.0 5.58 5.4 0.0 | |||
|
83 | g 9 1 10 10 10 10 10 10 10 10 10 -21 -7 9 6 10 10 10 10 10 9 9 9.4 9.76 9.4 0.0 | |||
|
84 | h 9 1 -20 -20 -20 -20 -20 -20 -20 -20 -21 -27 -24 -11 -18 -22 -21 -21 -21 -21 -20 -20 -19.7 -20.11 -20.5 0.0 | |||
|
85 | g 9 2 1 1 1 1 1 1 1 1 1 1 -1 -4 0 2 2 2 1 1 1 3 3.0 3.58 3.4 0.0 | |||
|
86 | h 9 2 14 14 14 14 14 14 14 15 15 17 19 12 12 15 16 16 16 15 15 15 13.4 12.69 11.6 0.0 | |||
|
87 | g 9 3 -11 -11 -11 -11 -11 -11 -12 -12 -12 -11 -25 -5 -9 -13 -12 -12 -12 -12 -12 -10 -8.4 -6.94 -5.3 0.0 | |||
|
88 | h 9 3 5 5 5 5 5 5 5 5 5 29 12 7 2 7 6 7 9 9 11 12 12.5 12.67 12.8 0.0 | |||
|
89 | g 9 4 12 12 12 12 12 12 12 11 11 3 10 2 1 10 10 10 9 9 9 8 6.3 5.01 3.1 0.0 | |||
|
90 | h 9 4 -3 -3 -3 -3 -3 -3 -3 -3 -3 -9 2 6 0 -4 -4 -4 -5 -6 -7 -6 -6.2 -6.72 -7.2 0.0 | |||
|
91 | g 9 5 1 1 1 1 1 1 1 1 1 16 5 4 4 -1 -1 -1 -3 -3 -4 -8 -8.9 -10.76 -12.4 0.0 | |||
|
92 | h 9 5 -2 -2 -2 -2 -2 -2 -2 -3 -3 4 2 -2 -3 -5 -5 -5 -6 -6 -7 -8 -8.4 -8.16 -7.4 0.0 | |||
|
93 | g 9 6 -2 -2 -2 -2 -2 -2 -2 -2 -2 -3 -5 1 -1 -1 0 -1 -1 -1 -2 -1 -1.5 -1.25 -0.8 0.0 | |||
|
94 | h 9 6 8 8 8 8 9 9 9 9 9 9 8 10 9 10 10 10 9 9 9 8 8.4 8.10 8.0 0.0 | |||
|
95 | g 9 7 2 2 2 2 2 2 3 3 3 -4 -2 2 -2 5 3 4 7 7 7 10 9.3 8.76 8.4 0.0 | |||
|
96 | h 9 7 10 10 10 10 10 10 10 11 11 6 8 7 8 10 11 11 10 9 8 5 3.8 2.92 2.2 0.0 | |||
|
97 | g 9 8 -1 0 0 0 0 0 0 0 1 -3 3 2 3 1 1 1 2 1 1 -2 -4.3 -6.66 -8.4 0.0 | |||
|
98 | h 9 8 -2 -2 -2 -2 -2 -2 -2 -2 -2 1 -11 -6 0 -4 -2 -3 -6 -7 -7 -8 -8.2 -7.73 -6.1 0.0 | |||
|
99 | g 9 9 -1 -1 -1 -1 -1 -1 -2 -2 -2 -4 8 5 -1 -2 -1 -2 -5 -5 -6 -8 -8.2 -9.22 -10.1 0.0 | |||
|
100 | h 9 9 2 2 2 2 2 2 2 2 2 8 -7 5 5 1 1 1 2 2 2 3 4.8 6.01 7.0 0.0 | |||
|
101 | g 10 0 -3 -3 -3 -3 -3 -3 -3 -3 -3 -3 -8 -3 1 -2 -3 -3 -4 -4 -3 -3 -2.6 -2.17 -2.0 0.0 | |||
|
102 | g 10 1 -4 -4 -4 -4 -4 -4 -4 -4 -4 11 4 -5 -3 -3 -3 -3 -4 -4 -4 -6 -6.0 -6.12 -6.3 0.0 | |||
|
103 | h 10 1 2 2 2 2 2 2 2 2 2 5 13 -4 4 2 1 1 1 1 2 1 1.7 2.19 2.8 0.0 | |||
|
104 | g 10 2 2 2 2 2 2 2 2 2 2 1 -1 -1 4 2 2 2 2 3 2 2 1.7 1.42 0.9 0.0 | |||
|
105 | h 10 2 1 1 1 1 1 1 1 1 1 1 -2 0 1 1 1 1 0 0 1 0 0.0 0.10 -0.1 0.0 | |||
|
106 | g 10 3 -5 -5 -5 -5 -5 -5 -5 -5 -5 2 13 2 0 -5 -5 -5 -5 -5 -5 -4 -3.1 -2.35 -1.1 0.0 | |||
|
107 | h 10 3 2 2 2 2 2 2 2 2 2 -20 -10 -8 0 2 3 3 3 3 3 4 4.0 4.46 4.7 0.0 | |||
|
108 | g 10 4 -2 -2 -2 -2 -2 -2 -2 -2 -2 -5 -4 -3 -1 -2 -1 -2 -2 -2 -2 -1 -0.5 -0.15 -0.2 0.0 | |||
|
109 | h 10 4 6 6 6 6 6 6 6 6 6 -1 2 -2 2 6 4 4 6 6 6 5 4.9 4.76 4.4 0.0 | |||
|
110 | g 10 5 6 6 6 6 6 6 6 6 6 -1 4 7 4 4 6 5 5 5 4 4 3.7 3.06 2.5 0.0 | |||
|
111 | h 10 5 -4 -4 -4 -4 -4 -4 -4 -4 -4 -6 -3 -4 -5 -4 -4 -4 -4 -4 -4 -5 -5.9 -6.58 -7.2 0.0 | |||
|
112 | g 10 6 4 4 4 4 4 4 4 4 4 8 12 4 6 4 4 4 3 3 3 2 1.0 0.29 -0.3 0.0 | |||
|
113 | h 10 6 0 0 0 0 0 0 0 0 0 6 6 1 1 0 0 -1 0 0 0 -1 -1.2 -1.01 -1.0 0.0 | |||
|
114 | g 10 7 0 0 0 0 0 0 0 0 0 -1 3 -2 1 0 1 1 1 1 1 2 2.0 2.06 2.2 0.0 | |||
|
115 | h 10 7 -2 -2 -2 -2 -2 -2 -2 -1 -1 -4 -3 -3 -1 -2 -1 -1 -1 -1 -2 -2 -2.9 -3.47 -4.0 0.0 | |||
|
116 | g 10 8 2 2 2 1 1 1 1 2 2 -3 2 6 -1 2 0 0 2 2 3 5 4.2 3.77 3.1 0.0 | |||
|
117 | h 10 8 4 4 4 4 4 4 4 4 4 -2 6 7 6 3 3 3 4 4 3 1 0.2 -0.86 -2.0 0.0 | |||
|
118 | g 10 9 2 2 2 2 3 3 3 3 3 5 10 -2 2 2 3 3 3 3 3 1 0.3 -0.21 -1.0 0.0 | |||
|
119 | h 10 9 0 0 0 0 0 0 0 0 0 0 11 -1 0 0 1 1 0 0 -1 -2 -2.2 -2.31 -2.0 0.0 | |||
|
120 | g 10 10 0 0 0 0 0 0 0 0 0 -2 3 0 0 0 -1 -1 0 0 0 0 -1.1 -2.09 -2.8 0.0 | |||
|
121 | h 10 10 -6 -6 -6 -6 -6 -6 -6 -6 -6 -2 8 -3 -7 -6 -4 -5 -6 -6 -6 -7 -7.4 -7.93 -8.3 0.0 | |||
|
122 | g 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2.7 2.95 3.0 0.0 | |||
|
123 | g 11 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1.7 -1.60 -1.5 0.0 | |||
|
124 | h 11 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.1 0.26 0.1 0.0 | |||
|
125 | g 11 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1.9 -1.88 -2.1 0.0 | |||
|
126 | h 11 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.3 1.44 1.7 0.0 | |||
|
127 | g 11 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.5 1.44 1.6 0.0 | |||
|
128 | h 11 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -0.77 -0.6 0.0 | |||
|
129 | g 11 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.1 -0.31 -0.5 0.0 | |||
|
130 | h 11 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2.6 -2.27 -1.8 0.0 | |||
|
131 | g 11 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.1 0.29 0.5 0.0 | |||
|
132 | h 11 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.9 0.90 0.9 0.0 | |||
|
133 | g 11 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.7 -0.79 -0.8 0.0 | |||
|
134 | h 11 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.7 -0.58 -0.4 0.0 | |||
|
135 | g 11 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.7 0.53 0.4 0.0 | |||
|
136 | h 11 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2.8 -2.69 -2.5 0.0 | |||
|
137 | g 11 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.7 1.80 1.8 0.0 | |||
|
138 | h 11 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -1.08 -1.3 0.0 | |||
|
139 | g 11 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.1 0.16 0.2 0.0 | |||
|
140 | h 11 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1.2 -1.58 -2.1 0.0 | |||
|
141 | g 11 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.2 0.96 0.8 0.0 | |||
|
142 | h 11 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1.9 -1.90 -1.9 0.0 | |||
|
143 | g 11 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4.0 3.99 3.8 0.0 | |||
|
144 | h 11 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -1.39 -1.8 0.0 | |||
|
145 | g 12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2.2 -2.15 -2.1 0.0 | |||
|
146 | g 12 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.3 -0.29 -0.2 0.0 | |||
|
147 | h 12 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.55 -0.8 0.0 | |||
|
148 | g 12 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.2 0.21 0.3 0.0 | |||
|
149 | h 12 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.23 0.3 0.0 | |||
|
150 | g 12 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.9 0.89 1.0 0.0 | |||
|
151 | h 12 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2.5 2.38 2.2 0.0 | |||
|
152 | g 12 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.2 -0.38 -0.7 0.0 | |||
|
153 | h 12 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2.6 -2.63 -2.5 0.0 | |||
|
154 | g 12 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.9 0.96 0.9 0.0 | |||
|
155 | h 12 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.7 0.61 0.5 0.0 | |||
|
156 | g 12 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.5 -0.30 -0.1 0.0 | |||
|
157 | h 12 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.40 0.6 0.0 | |||
|
158 | g 12 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.46 0.5 0.0 | |||
|
159 | h 12 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.0 0.01 0.0 0.0 | |||
|
160 | g 12 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.3 -0.35 -0.4 0.0 | |||
|
161 | h 12 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.0 0.02 0.1 0.0 | |||
|
162 | g 12 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.36 -0.4 0.0 | |||
|
163 | h 12 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.28 0.3 0.0 | |||
|
164 | g 12 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.1 0.08 0.2 0.0 | |||
|
165 | h 12 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -0.87 -0.9 0.0 | |||
|
166 | g 12 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.2 -0.49 -0.8 0.0 | |||
|
167 | h 12 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.34 -0.2 0.0 | |||
|
168 | g 12 12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.08 0.0 0.0 | |||
|
169 | h 12 12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.8 0.88 0.8 0.0 | |||
|
170 | g 13 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.2 -0.16 -0.2 0.0 | |||
|
171 | g 13 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -0.88 -0.9 0.0 | |||
|
172 | h 13 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -0.76 -0.8 0.0 | |||
|
173 | g 13 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.30 0.3 0.0 | |||
|
174 | h 13 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.2 0.33 0.3 0.0 | |||
|
175 | g 13 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.1 0.28 0.4 0.0 | |||
|
176 | h 13 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.8 1.72 1.7 0.0 | |||
|
177 | g 13 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.43 -0.4 0.0 | |||
|
178 | h 13 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.54 -0.6 0.0 | |||
|
179 | g 13 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.3 1.18 1.1 0.0 | |||
|
180 | h 13 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1.0 -1.07 -1.2 0.0 | |||
|
181 | g 13 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.37 -0.3 0.0 | |||
|
182 | h 13 6 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.1 -0.04 -0.1 0.0 | |||
|
183 | g 13 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.7 0.75 0.8 0.0 | |||
|
184 | h 13 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.7 0.63 0.5 0.0 | |||
|
185 | g 13 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.4 -0.26 -0.2 0.0 | |||
|
186 | h 13 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.21 0.1 0.0 | |||
|
187 | g 13 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.35 0.4 0.0 | |||
|
188 | h 13 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.6 0.53 0.5 0.0 | |||
|
189 | g 13 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.1 -0.05 0.0 0.0 | |||
|
190 | h 13 10 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.3 0.38 0.4 0.0 | |||
|
191 | g 13 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.4 0.41 0.4 0.0 | |||
|
192 | h 13 11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.2 -0.22 -0.2 0.0 | |||
|
193 | g 13 12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.0 -0.10 -0.3 0.0 | |||
|
194 | h 13 12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.5 -0.57 -0.5 0.0 | |||
|
195 | g 13 13 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.1 -0.18 -0.3 0.0 | |||
|
196 | h 13 13 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0.9 -0.82 -0.8 0.0 |
1 | NO CONTENT: new file 100644 |
|
NO CONTENT: new file 100644 | ||
The requested commit or file is too big and content was truncated. Show full diff |
1 | NO CONTENT: new file 100644, binary diff hidden |
|
NO CONTENT: new file 100644, binary diff hidden |
@@ -148,9 +148,3 h2 { | |||||
148 | .container, .container-lg, .container-md, .container-sm, .container-xl { |
|
148 | .container, .container-lg, .container-md, .container-sm, .container-xl { | |
149 | max-width: 1020px; |
|
149 | max-width: 1020px; | |
150 | } |
|
150 | } | |
151 |
|
||||
152 | @media (min-width: 1380px) { |
|
|||
153 | .card-columns { |
|
|||
154 | column-count: 5; |
|
|||
155 | } |
|
|||
156 | } |
|
@@ -46,57 +46,58 | |||||
46 | <div class="invalid-tooltip"> |
|
46 | <div class="invalid-tooltip"> | |
47 | Please enter a valid date. |
|
47 | Please enter a valid date. | |
48 | </div> |
|
48 | </div> | |
49 |
<select |
|
49 | <select id="overjro-experiment" class="form-control form-control-sm"> | |
50 | <option value="-1">Experiment:</option> |
|
50 | <option value="-1">Experiment:</option> | |
51 | <option value="-1">------------------</option> |
|
51 | <option value="-1">------------------</option> | |
52 |
<option value=" |
|
52 | <option value="50">Vertical Drifts</option> | |
53 |
<option value=" |
|
53 | <option value="51">East West 1996 (W beam)</option> | |
54 |
<option value=" |
|
54 | <option value="52">East West 1996 (E beam)</option> | |
55 |
<option value=" |
|
55 | <option value="61">East West 2003</option> | |
56 |
<option value=" |
|
56 | <option value="60">Differential Phase 2000</option> | |
57 |
<option value=" |
|
57 | <option value="63">Differential Phase 2004 High Alt</option> | |
58 |
<option value=" |
|
58 | <option value="64">Differential Phase 2005 - 2006</option> | |
59 |
|
|
59 | <option value="54">DEWD 2005</option> | |
|
60 | <option value="53">DVD 2006 - 2008</option> | |||
60 | <option value="-1">------------------</option> |
|
61 | <option value="-1">------------------</option> | |
61 |
<option value=" |
|
62 | <option value="4">Oblique ISR On-Axis</option> | |
62 |
<option value=" |
|
63 | <option value="5">Oblique ISR 4.5</option> | |
63 |
<option value=" |
|
64 | <option value="6">Oblique ISR 6.0S</option> | |
64 |
<option value=" |
|
65 | <option value="7">Oblique ISR 3.0N</option> | |
65 | <option value="-1">------------------</option> |
|
66 | <option value="-1">------------------</option> | |
66 |
<option value=" |
|
67 | <option value="16">JULIA CP2</option> | |
67 |
<option value=" |
|
68 | <option value="17">JULIA CP3</option> | |
68 |
<option value=" |
|
69 | <option value="18">JULIA V (2005-2006)</option> | |
69 |
<option value=" |
|
70 | <option value="65">JULIA EW 2003</option> | |
70 |
|
|
71 | <option value="19">JULIA EW (2006-2007)</option> | |
71 | <option value="-1">------------------</option> |
|
72 | <option value="-1">------------------</option> | |
72 | <option value="0">Modulo Rx</option> |
|
73 | <option value="0">Modulo Rx</option> | |
73 | <option value="1">1/16 Rx</option> |
|
74 | <option value="1">1/16 Rx</option> | |
74 | <option value="2">1/4 Rx</option> |
|
75 | <option value="2">1/4 Rx</option> | |
75 | <option value="3">All Rx</option> |
|
76 | <option value="3">All Rx</option> | |
76 | <option value="-1">------------------</option> |
|
77 | <option value="-1">------------------</option> | |
77 |
<option value=" |
|
78 | <option value="21">EW Imaging 1996</option> | |
78 |
<option value=" |
|
79 | <option value="22">EW Imaging 2003</option> | |
79 |
|
|
80 | <option value="23">EW Imaging 2006-2008</option> | |
80 | <option value="-1">------------------</option> |
|
81 | <option value="-1">------------------</option> | |
81 |
<option value=" |
|
82 | <option value="359">MST North (Fritts)</option> | |
82 |
<option value=" |
|
83 | <option value="360">MST West (Fritts)</option> | |
83 |
<option value=" |
|
84 | <option value="361">MST South (Fritts)</option> | |
84 |
<option value=" |
|
85 | <option value="362">MST East (Fritts)</option> | |
85 | <option value="-1">------------------</option> |
|
86 | <option value="-1">------------------</option> | |
86 |
<option value=" |
|
87 | <option value="67">Vertical (Yellow Cables)</option> | |
87 | </select> |
|
88 | </select> | |
88 | <br> |
|
89 | <br> | |
89 | <p class="card-text">Choose object: |
|
90 | <p class="card-text">Choose object: | |
90 | <div class="form-check card-text"> |
|
91 | <div class="form-check card-text"> | |
91 |
<input class="form-check-input" type="checkbox" id="inlineCheckbox1" value=" |
|
92 | <input class="form-check-input" type="checkbox" id="inlineCheckbox1" value="bfield" name="celestial"> | |
92 | <label class="form-check-label" for="inlineCheckbox1">B Field</label><br> |
|
93 | <label class="form-check-label" for="inlineCheckbox1">B Field</label><br> | |
93 |
<input class="form-check-input" type="checkbox" id="inlineCheckbox1" value=" |
|
94 | <input class="form-check-input" type="checkbox" id="inlineCheckbox1" value="sun" name="celestial"> | |
94 | <label class="form-check-label" for="inlineCheckbox1">Sun</label><br> |
|
95 | <label class="form-check-label" for="inlineCheckbox1">Sun</label><br> | |
95 |
<input class="form-check-input" type="checkbox" id="inlineCheckbox1" value=" |
|
96 | <input class="form-check-input" type="checkbox" id="inlineCheckbox1" value="moon" name="celestial"> | |
96 | <label class="form-check-label" for="inlineCheckbox1">Moon</label><br> |
|
97 | <label class="form-check-label" for="inlineCheckbox1">Moon</label><br> | |
97 |
<input class="form-check-input" type="checkbox" id="inlineCheckbox1" value=" |
|
98 | <input class="form-check-input" type="checkbox" id="inlineCheckbox1" value="hydra" name="celestial"> | |
98 | <label class="form-check-label" for="inlineCheckbox1">Hydra</label><br> |
|
99 | <label class="form-check-label" for="inlineCheckbox1">Hydra</label><br> | |
99 |
<input class="form-check-input" type="checkbox" id="inlineCheckbox1" value=" |
|
100 | <input class="form-check-input" type="checkbox" id="inlineCheckbox1" value="galaxy" name="celestial"> | |
100 | <label class="form-check-label" for="inlineCheckbox1">Galaxy Center</label> |
|
101 | <label class="form-check-label" for="inlineCheckbox1">Galaxy Center</label> | |
101 | </div> |
|
102 | </div> | |
102 | <br> |
|
103 | <br> | |
@@ -106,9 +107,9 | |||||
106 | value="5.0" required> |
|
107 | value="5.0" required> | |
107 | </div> |
|
108 | </div> | |
108 | <div class="form-group card-text"> |
|
109 | <div class="form-group card-text"> | |
109 | <label class="form-check-label" for="overjro-height">Height [km]:</label> |
|
110 | <label class="form-check-label" for="overjro-height">Heights [km]:</label> | |
110 | <input type="text" class="form-control form-control-sm tools-date" id="overjro-height" placeholder="Enter Height" |
|
111 | <input type="text" class="form-control form-control-sm tools-date" id="overjro-height" placeholder="Enter Heights [km]" | |
111 | value="" required> |
|
112 | value="100" required> | |
112 | </div> |
|
113 | </div> | |
113 | </p> |
|
114 | </p> | |
114 | </div> |
|
115 | </div> | |
@@ -142,6 +143,9 | |||||
142 | <div class="modal-body text-center"> |
|
143 | <div class="modal-body text-center"> | |
143 | <img class="img-fluid" src=""> |
|
144 | <img class="img-fluid" src=""> | |
144 | </div> |
|
145 | </div> | |
|
146 | <div class="modal-body text-center"> | |||
|
147 | <p></p> | |||
|
148 | </div> | |||
145 | </div> |
|
149 | </div> | |
146 | </div> |
|
150 | </div> | |
147 | </div> |
|
151 | </div> | |
@@ -156,15 +160,39 | |||||
156 | //get data attribute of the clicked element |
|
160 | //get data attribute of the clicked element | |
157 | var title = $(e.relatedTarget).data('title'); |
|
161 | var title = $(e.relatedTarget).data('title'); | |
158 | var image = $(e.relatedTarget).data('image'); |
|
162 | var image = $(e.relatedTarget).data('image'); | |
|
163 | $(e.currentTarget).find('p').text(''); | |||
|
164 | $(e.currentTarget).find('img').attr('src', ''); | |||
159 |
|
165 | |||
160 | if (image.indexOf('skynoise') > 0) { |
|
166 | if (image.indexOf('skynoise') > 0) { | |
161 | var dt = $('#skynoise-date').val(); |
|
167 | var dt = $('#skynoise-date').val(); | |
162 | image += '?date=' + dt; |
|
168 | image += '?date=' + dt; | |
|
169 | //populate values | |||
|
170 | $(e.currentTarget).find('h5').text(title); | |||
|
171 | $(e.currentTarget).find('img').attr('src', image); | |||
163 | } |
|
172 | } | |
164 |
|
173 | |||
165 | //populate values |
|
174 | if (image.indexOf('overjro') > 0) { | |
166 | $(e.currentTarget).find('h5').text(title); |
|
175 | $(e.currentTarget).find('h5').text(title); | |
|
176 | ||||
|
177 | if ($('#overjro-experiment').val() == '-1'){ | |||
|
178 | $(e.currentTarget).find('p').text('Missing Experiment'); | |||
|
179 | } else { | |||
|
180 | ||||
|
181 | var dt = $('#overjro-date').val(); | |||
|
182 | var favorite = []; | |||
|
183 | $.each($("input[name='celestial']:checked"), function(){ | |||
|
184 | favorite.push($(this).val()); | |||
|
185 | }); | |||
|
186 | ||||
|
187 | image += '?date=' + dt; | |||
|
188 | image += '&experiment=' + $('#overjro-experiment').val(); | |||
|
189 | image += '&angle=' + $('#overjro-angle').val(); | |||
|
190 | image += '&height=' + $('#overjro-height').val(); | |||
|
191 | image += '&bodys=' + favorite.join(","); | |||
|
192 | ||||
167 | $(e.currentTarget).find('img').attr('src', image); |
|
193 | $(e.currentTarget).find('img').attr('src', image); | |
|
194 | } | |||
|
195 | } | |||
168 | }); |
|
196 | }); | |
169 |
|
197 | |||
170 | $('#doy-date').change(function() { |
|
198 | $('#doy-date').change(function() { |
@@ -16,7 +16,7 import mongoengine | |||||
16 |
|
16 | |||
17 | from plotter.models import Experiment, ExpDetail, PlotMeta, PlotData, JROReport |
|
17 | from plotter.models import Experiment, ExpDetail, PlotMeta, PlotData, JROReport | |
18 |
|
18 | |||
19 | from utils.plots import skynoise_plot |
|
19 | from utils.plots import skynoise_plot, overjro_plot | |
20 |
|
20 | |||
21 | host = os.environ.get('HOST_MONGO', 'localhost') |
|
21 | host = os.environ.get('HOST_MONGO', 'localhost') | |
22 | mongoengine.connect('dbplots', host=host, port=27017) |
|
22 | mongoengine.connect('dbplots', host=host, port=27017) | |
@@ -258,7 +258,12 def plot_overjro(request): | |||||
258 | else: |
|
258 | else: | |
259 | date = datetime.strptime(date, '%Y-%m-%d') |
|
259 | date = datetime.strptime(date, '%Y-%m-%d') | |
260 |
|
260 | |||
261 | data = skynoise_plot(date.year, date.month, date.day) |
|
261 | pattern = int(request.GET.get('experiment', '1')) | |
|
262 | angle = float(request.GET.get('angle', '5')) | |||
|
263 | height = [float(h) for h in request.GET.get('height', '100').split(',')] | |||
|
264 | bodys = (request.GET.get('bodys', '')).split(',') | |||
|
265 | ||||
|
266 | data = overjro_plot(pattern, date, angle, height, bodys) | |||
262 | response = HttpResponse(data.getvalue(), content_type='image/png') |
|
267 | response = HttpResponse(data.getvalue(), content_type='image/png') | |
263 |
|
268 | |||
264 | return response |
|
269 | return response |
This diff has been collapsed as it changes many lines, (1661 lines changed) Show them Hide them | |||||
@@ -4,11 +4,1637 import time | |||||
4 | import numpy |
|
4 | import numpy | |
5 | import scipy |
|
5 | import scipy | |
6 | import base64 |
|
6 | import base64 | |
|
7 | import datetime | |||
7 | from io import BytesIO |
|
8 | from io import BytesIO | |
8 | import scipy.interpolate |
|
9 | import scipy.interpolate | |
9 | from matplotlib.figure import Figure |
|
10 | from matplotlib.figure import Figure | |
10 |
|
11 | |||
11 | from utils import TimeTools |
|
12 | from utils import TimeTools | |
|
13 | from utils.patterns import select_pattern | |||
|
14 | from utils import Misc_Routines | |||
|
15 | from utils import Astro_Coords | |||
|
16 | from utils import gaussfit | |||
|
17 | ||||
|
18 | ||||
|
19 | attenuation = numpy.array([[[-21.25,-15.25,-9.25,-3.25,3.25,9.25,15.25,21.25], | |||
|
20 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
21 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
22 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
23 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
24 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
25 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25], | |||
|
26 | [-21.25,-15.25,-9.25,-3.25,03.25,09.25,15.25,21.25]], | |||
|
27 | [[21.25,21.25,21.25,21.25,21.25,21.25,21.25,21.25], | |||
|
28 | [15.25,15.25,15.25,15.25,15.25,15.25,15.25,15.25], | |||
|
29 | [09.25,09.25,09.25,09.25,09.25,09.25,09.25,09.25], | |||
|
30 | [03.25,03.25,03.25,03.25,03.25,03.25,03.25,03.25], | |||
|
31 | [-03.25,-03.25,-03.25,-03.25,-03.25,-03.25,-03.25,-03.25], | |||
|
32 | [-09.25,-09.25,-09.25,-09.25,-09.25,-09.25,-09.25,-09.25], | |||
|
33 | [-15.25,-15.25,-15.25,-15.25,-15.25,-15.25,-15.25,-15.25], | |||
|
34 | [-21.25,-21.25,-21.25,-21.25,-21.25,-21.25,-21.25,-21.25]]]) | |||
|
35 | ||||
|
36 | class BField(): | |||
|
37 | def __init__(self,year=None,doy=None,site=1,heights=None,alpha_i=90): | |||
|
38 | """ | |||
|
39 | BField class creates an object to get the Magnetic field for a specific date and | |||
|
40 | height(s). | |||
|
41 | ||||
|
42 | Parameters | |||
|
43 | ---------- | |||
|
44 | year = A scalar giving the desired year. If the value is None (default value) then | |||
|
45 | the current year will be used. | |||
|
46 | doy = A scalar giving the desired day of the year. If the value is None (default va- | |||
|
47 | lue) then the current doy will be used. | |||
|
48 | site = An integer to choose the geographic coordinates of the place where the magne- | |||
|
49 | tic field will be computed. The default value is over Jicamarca (site=1) | |||
|
50 | heights = An array giving the heights (km) where the magnetic field will be modeled By default the magnetic field will be computed at 100, 500 and 1000km. | |||
|
51 | alpha_i = Angle to interpolate the magnetic field. | |||
|
52 | ||||
|
53 | Modification History | |||
|
54 | -------------------- | |||
|
55 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 07 October 2009. | |||
|
56 | """ | |||
|
57 | ||||
|
58 | tmp = time.localtime() | |||
|
59 | if year==None: year = tmp[0] | |||
|
60 | if doy==None: doy = tmp[7] | |||
|
61 | self.year = year | |||
|
62 | self.doy = doy | |||
|
63 | self.site = site | |||
|
64 | if heights is None: | |||
|
65 | heights = numpy.array([100,500,1000]) | |||
|
66 | else: | |||
|
67 | heights = numpy.array(heights) | |||
|
68 | self.heights = heights | |||
|
69 | self.alpha_i = alpha_i | |||
|
70 | ||||
|
71 | def getBField(self,maglimits=numpy.array([-7,-7,7,7])): | |||
|
72 | """ | |||
|
73 | getBField models the magnetic field for a different heights in a specific date. | |||
|
74 | ||||
|
75 | Parameters | |||
|
76 | ---------- | |||
|
77 | maglimits = An 4-elements array giving ..... The default value is [-7,-7,7,7]. | |||
|
78 | ||||
|
79 | Return | |||
|
80 | ------ | |||
|
81 | dcos = An 4-dimensional array giving the directional cosines of the magnetic field | |||
|
82 | over the desired place. | |||
|
83 | alpha = An 3-dimensional array giving the angle of the magnetic field over the desi- | |||
|
84 | red place. | |||
|
85 | ||||
|
86 | Modification History | |||
|
87 | -------------------- | |||
|
88 | Converted to Python by Freddy R. Galindo, ROJ, 07 October 2009. | |||
|
89 | """ | |||
|
90 | ||||
|
91 | x_ant = numpy.array([1,0,0]) | |||
|
92 | y_ant = numpy.array([0,1,0]) | |||
|
93 | z_ant = numpy.array([0,0,1]) | |||
|
94 | ||||
|
95 | if self.site==0: | |||
|
96 | title_site = "Magnetic equator" | |||
|
97 | coord_site = numpy.array([-76+52./60.,-11+57/60.,0.5]) | |||
|
98 | elif self.site==1: | |||
|
99 | title_site = 'Jicamarca' | |||
|
100 | coord_site = [-76-52./60.,-11-57/60.,0.5] | |||
|
101 | theta = (45+5.35)*numpy.pi/180. # (50.35 and 1.46 from Fleish Thesis) | |||
|
102 | delta = -1.46*numpy.pi/180 | |||
|
103 | ||||
|
104 | x_ant1 = numpy.roll(self.rotvector(self.rotvector(x_ant,1,delta),3,theta),1) | |||
|
105 | y_ant1 = numpy.roll(self.rotvector(self.rotvector(y_ant,1,delta),3,theta),1) | |||
|
106 | z_ant1 = numpy.roll(self.rotvector(self.rotvector(z_ant,1,delta),3,theta),1) | |||
|
107 | ||||
|
108 | ang0 = -1*coord_site[0]*numpy.pi/180. | |||
|
109 | ang1 = coord_site[1]*numpy.pi/180. | |||
|
110 | x_ant = self.rotvector(self.rotvector(x_ant1,2,ang1),3,ang0) | |||
|
111 | y_ant = self.rotvector(self.rotvector(y_ant1,2,ang1),3,ang0) | |||
|
112 | z_ant = self.rotvector(self.rotvector(z_ant1,2,ang1),3,ang0) | |||
|
113 | else: | |||
|
114 | # print "No defined Site. Skip..." | |||
|
115 | return None | |||
|
116 | ||||
|
117 | nhei = self.heights.size | |||
|
118 | pt_intercep = numpy.zeros((nhei,2)) | |||
|
119 | nfields = 1 | |||
|
120 | ||||
|
121 | grid_res = 0.5 | |||
|
122 | nlon = int(numpy.int(maglimits[2] - maglimits[0])/grid_res + 1) | |||
|
123 | nlat = int(numpy.int(maglimits[3] - maglimits[1])/grid_res + 1) | |||
|
124 | ||||
|
125 | location = numpy.zeros((nlon,nlat,2)) | |||
|
126 | mlon = numpy.atleast_2d(numpy.arange(nlon)*grid_res + maglimits[0]) | |||
|
127 | mrep = numpy.atleast_2d(numpy.zeros(nlat) + 1) | |||
|
128 | location0 = numpy.dot(mlon.transpose(),mrep) | |||
|
129 | ||||
|
130 | mlat = numpy.atleast_2d(numpy.arange(nlat)*grid_res + maglimits[1]) | |||
|
131 | mrep = numpy.atleast_2d(numpy.zeros(nlon) + 1) | |||
|
132 | location1 = numpy.dot(mrep.transpose(),mlat) | |||
|
133 | ||||
|
134 | location[:,:,0] = location0 | |||
|
135 | location[:,:,1] = location1 | |||
|
136 | ||||
|
137 | alpha = numpy.zeros((nlon,nlat,nhei)) | |||
|
138 | rr = numpy.zeros((nlon,nlat,nhei,3)) | |||
|
139 | dcos = numpy.zeros((nlon,nlat,nhei,2)) | |||
|
140 | ||||
|
141 | global first_time | |||
|
142 | ||||
|
143 | first_time = None | |||
|
144 | for ilon in numpy.arange(nlon): | |||
|
145 | for ilat in numpy.arange(nlat): | |||
|
146 | outs = self.__bdotk(self.heights, | |||
|
147 | self.year + self.doy/366., | |||
|
148 | coord_site[1], | |||
|
149 | coord_site[0], | |||
|
150 | coord_site[2], | |||
|
151 | coord_site[1]+location[ilon,ilat,1], | |||
|
152 | location[ilon,ilat,0]*720./180.) | |||
|
153 | ||||
|
154 | alpha[ilon, ilat,:] = outs[1] | |||
|
155 | rr[ilon, ilat,:,:] = outs[3] | |||
|
156 | ||||
|
157 | mrep = numpy.atleast_2d((numpy.zeros(nhei)+1)).transpose() | |||
|
158 | tmp = outs[3]*numpy.dot(mrep,numpy.atleast_2d(x_ant)) | |||
|
159 | tmp = tmp.sum(axis=1) | |||
|
160 | dcos[ilon,ilat,:,0] = tmp/numpy.sqrt((outs[3]**2).sum(axis=1)) | |||
|
161 | ||||
|
162 | mrep = numpy.atleast_2d((numpy.zeros(nhei)+1)).transpose() | |||
|
163 | tmp = outs[3]*numpy.dot(mrep,numpy.atleast_2d(y_ant)) | |||
|
164 | tmp = tmp.sum(axis=1) | |||
|
165 | dcos[ilon,ilat,:,1] = tmp/numpy.sqrt((outs[3]**2).sum(axis=1)) | |||
|
166 | ||||
|
167 | return dcos, alpha, nlon, nlat | |||
|
168 | ||||
|
169 | ||||
|
170 | def __bdotk(self,heights,tm,gdlat=-11.95,gdlon=-76.8667,gdalt=0.0,decd=-12.88, ham=-4.61666667): | |||
|
171 | ||||
|
172 | global first_time | |||
|
173 | # Mean Earth radius in Km WGS 84 | |||
|
174 | a_igrf = 6371.2 | |||
|
175 | ||||
|
176 | bk = numpy.zeros(heights.size) | |||
|
177 | alpha = numpy.zeros(heights.size) | |||
|
178 | bfm = numpy.zeros(heights.size) | |||
|
179 | rr = numpy.zeros((heights.size,3)) | |||
|
180 | rgc = numpy.zeros((heights.size,3)) | |||
|
181 | ||||
|
182 | ObjGeodetic = Astro_Coords.Geodetic(gdlat,gdalt) | |||
|
183 | [gclat,gcalt] = ObjGeodetic.change2geocentric() | |||
|
184 | ||||
|
185 | gclat = gclat*numpy.pi/180. | |||
|
186 | gclon = gdlon*numpy.pi/180. | |||
|
187 | ||||
|
188 | # Antenna position from center of Earth | |||
|
189 | ca_vector = [numpy.cos(gclat)*numpy.cos(gclon),numpy.cos(gclat)*numpy.sin(gclon),numpy.sin(gclat)] | |||
|
190 | ca_vector = gcalt*numpy.array(ca_vector) | |||
|
191 | ||||
|
192 | dec = decd*numpy.pi/180. | |||
|
193 | ||||
|
194 | # K vector respect to the center of earth. | |||
|
195 | klon = gclon + ham*numpy.pi/720. | |||
|
196 | k_vector = [numpy.cos(dec)*numpy.cos(klon),numpy.cos(dec)*numpy.sin(klon),numpy.sin(dec)] | |||
|
197 | k_vector = numpy.array(k_vector) | |||
|
198 | ||||
|
199 | for ih in numpy.arange(heights.size): | |||
|
200 | # Vector from Earth's center to volume of interest | |||
|
201 | rr[ih,:] = k_vector*heights[ih] | |||
|
202 | cv_vector = numpy.squeeze(ca_vector) + rr[ih,:] | |||
|
203 | ||||
|
204 | cv_gcalt = numpy.sqrt(numpy.sum(cv_vector**2.)) | |||
|
205 | cvxy = numpy.sqrt(numpy.sum(cv_vector[0:2]**2.)) | |||
|
206 | ||||
|
207 | radial = cv_vector/cv_gcalt | |||
|
208 | east = numpy.array([-1*cv_vector[1],cv_vector[0],0])/cvxy | |||
|
209 | comp1 = east[1]*radial[2] - radial[1]*east[2] | |||
|
210 | comp2 = east[2]*radial[0] - radial[2]*east[0] | |||
|
211 | comp3 = east[0]*radial[1] - radial[0]*east[1] | |||
|
212 | north = -1*numpy.array([comp1, comp2, comp3]) | |||
|
213 | ||||
|
214 | rr_k = cv_vector - numpy.squeeze(ca_vector) | |||
|
215 | u_rr = rr_k/numpy.sqrt(numpy.sum(rr_k**2.)) | |||
|
216 | ||||
|
217 | cv_gclat = numpy.arctan2(cv_vector[2],cvxy) | |||
|
218 | cv_gclon = numpy.arctan2(cv_vector[1],cv_vector[0]) | |||
|
219 | ||||
|
220 | bhei = cv_gcalt-a_igrf | |||
|
221 | blat = cv_gclat*180./numpy.pi | |||
|
222 | blon = cv_gclon*180./numpy.pi | |||
|
223 | bfield = self.__igrfkudeki(bhei,tm,blat,blon) | |||
|
224 | ||||
|
225 | B = (bfield[0]*north + bfield[1]*east - bfield[2]*radial)*1.0e-5 | |||
|
226 | ||||
|
227 | bfm[ih] = numpy.sqrt(numpy.sum(B**2.)) #module | |||
|
228 | bk[ih] = numpy.sum(u_rr*B) | |||
|
229 | alpha[ih] = numpy.arccos(bk[ih]/bfm[ih])*180/numpy.pi | |||
|
230 | rgc[ih,:] = numpy.array([cv_gclon, cv_gclat, cv_gcalt]) | |||
|
231 | ||||
|
232 | return bk, alpha, bfm, rr, rgc | |||
|
233 | ||||
|
234 | ||||
|
235 | def __igrfkudeki(self,heights,time,latitude,longitude,ae=6371.2): | |||
|
236 | """ | |||
|
237 | __igrfkudeki calculates the International Geomagnetic Reference Field for given in- | |||
|
238 | put conditions based on IGRF2005 coefficients. | |||
|
239 | ||||
|
240 | Parameters | |||
|
241 | ---------- | |||
|
242 | heights = Scalar or vector giving the height above the Earth of the point in ques- | |||
|
243 | tion in kilometers. | |||
|
244 | time = Scalar or vector giving the decimal year of time in question (e.g. 1991.2). | |||
|
245 | latitude = Latitude of point in question in decimal degrees. Scalar or vector. | |||
|
246 | longitude = Longitude of point in question in decimal degrees. Scalar or vector. | |||
|
247 | ae = | |||
|
248 | first_time = | |||
|
249 | ||||
|
250 | Return | |||
|
251 | ------ | |||
|
252 | bn = | |||
|
253 | be = | |||
|
254 | bd = | |||
|
255 | bmod = | |||
|
256 | balpha = | |||
|
257 | first_time = | |||
|
258 | ||||
|
259 | Modification History | |||
|
260 | -------------------- | |||
|
261 | Converted to Python by Freddy R. Galindo, ROJ, 03 October 2009. | |||
|
262 | """ | |||
|
263 | ||||
|
264 | global first_time | |||
|
265 | global gs, hs, nvec, mvec, maxcoef | |||
|
266 | ||||
|
267 | heights = numpy.atleast_1d(heights) | |||
|
268 | time = numpy.atleast_1d(time) | |||
|
269 | latitude = numpy.atleast_1d(latitude) | |||
|
270 | longitude = numpy.atleast_1d(longitude) | |||
|
271 | ||||
|
272 | if numpy.max(latitude)==90: | |||
|
273 | # print "Field calculations are not supported at geographic poles" | |||
|
274 | pass | |||
|
275 | ||||
|
276 | # output arrays | |||
|
277 | bn = numpy.zeros(heights.size) | |||
|
278 | be = numpy.zeros(heights.size) | |||
|
279 | bd = numpy.zeros(heights.size) | |||
|
280 | ||||
|
281 | if first_time==None:first_time=0 | |||
|
282 | ||||
|
283 | time0 = time[0] | |||
|
284 | if time!=first_time: | |||
|
285 | #print "Getting coefficients for", time0 | |||
|
286 | [periods,g,h ] = self.__readIGRFcoeff() | |||
|
287 | top_year = numpy.max(periods) | |||
|
288 | nperiod = (top_year - 1900)/5 + 1 | |||
|
289 | ||||
|
290 | maxcoef = 10 | |||
|
291 | if time0>=2000:maxcoef = 12 | |||
|
292 | ||||
|
293 | ||||
|
294 | # Normalization array for Schmidt fucntions | |||
|
295 | multer = numpy.zeros((2+maxcoef,1+maxcoef)) + 1 | |||
|
296 | for cn in (numpy.arange(maxcoef)+1): | |||
|
297 | for rm in (numpy.arange(cn)+1): | |||
|
298 | tmp = numpy.arange(2*rm) + cn - rm + 1. | |||
|
299 | multer[rm+1,cn] = ((-1.)**rm)*numpy.sqrt(2./tmp.prod()) | |||
|
300 | ||||
|
301 | schmidt = multer[1:,1:].transpose() | |||
|
302 | ||||
|
303 | # n and m arrays | |||
|
304 | nvec = numpy.atleast_2d(numpy.arange(maxcoef)+2) | |||
|
305 | mvec = numpy.atleast_2d(numpy.arange(maxcoef+1)).transpose() | |||
|
306 | ||||
|
307 | # Time adjusted igrf g and h with Schmidt normalization | |||
|
308 | # IGRF coefficient arrays: g0(n,m), n=1, maxcoeff,m=0, maxcoeff, ... | |||
|
309 | if time0<top_year: | |||
|
310 | dtime = (time0 - 1900) % 5 | |||
|
311 | ntime = (time0 - 1900 - dtime)/5 | |||
|
312 | else: | |||
|
313 | # Estimating coefficients for times > top_year | |||
|
314 | dtime = (time0 - top_year) + 5 | |||
|
315 | ntime = g[:,0,0].size - 2 | |||
|
316 | ||||
|
317 | g0 = g[ntime,1:maxcoef+1,:maxcoef+1] | |||
|
318 | h0 = h[ntime,1:maxcoef+1,:maxcoef+1] | |||
|
319 | gdot = g[ntime+1,1:maxcoef+1,:maxcoef+1]-g[ntime,1:maxcoef+1,:maxcoef+1] | |||
|
320 | hdot = h[ntime+1,1:maxcoef+1,:maxcoef+1]-h[ntime,1:maxcoef+1,:maxcoef+1] | |||
|
321 | gs = (g0 + dtime*(gdot/5.))*schmidt[:maxcoef,0:maxcoef+1] | |||
|
322 | hs = (h0 + dtime*(hdot/5.))*schmidt[:maxcoef,0:maxcoef+1] | |||
|
323 | ||||
|
324 | first_time = time0 | |||
|
325 | ||||
|
326 | for ii in numpy.arange(heights.size): | |||
|
327 | # Height dependence array rad = (ae/(ae+height))**(n+3) | |||
|
328 | rad = numpy.atleast_2d((ae/(ae + heights[ii]))**(nvec+1)) | |||
|
329 | ||||
|
330 | # Sin and Cos of m times longitude phi arrays | |||
|
331 | mphi = mvec*longitude[ii]*numpy.pi/180. | |||
|
332 | cosmphi = numpy.atleast_2d(numpy.cos(mphi)) | |||
|
333 | sinmphi = numpy.atleast_2d(numpy.sin(mphi)) | |||
|
334 | ||||
|
335 | # Cos of colatitude theta | |||
|
336 | c = numpy.cos((90 - latitude[ii])*numpy.pi/180.) | |||
|
337 | ||||
|
338 | # Legendre functions p(n,m|c) | |||
|
339 | [p,dp]= scipy.special.lpmn(maxcoef+1,maxcoef+1,c) | |||
|
340 | p = p[:,:-1].transpose() | |||
|
341 | s = numpy.sqrt((1. - c)*(1 + c)) | |||
|
342 | ||||
|
343 | # Generate derivative array dpdtheta = -s*dpdc | |||
|
344 | dpdtheta = c*p/s | |||
|
345 | for m in numpy.arange(maxcoef+2): dpdtheta[:,m] = m*dpdtheta[:,m] | |||
|
346 | dpdtheta = dpdtheta + numpy.roll(p,-1,axis=1) | |||
|
347 | ||||
|
348 | # Extracting arrays required for field calculations | |||
|
349 | p = p[1:maxcoef+1,:maxcoef+1] | |||
|
350 | dpdtheta = dpdtheta[1:maxcoef+1,:maxcoef+1] | |||
|
351 | ||||
|
352 | # Weigh p and dpdtheta with gs and hs coefficients. | |||
|
353 | gp = gs*p | |||
|
354 | hp = hs*p | |||
|
355 | gdpdtheta = gs*dpdtheta | |||
|
356 | hdpdtheta = hs*dpdtheta | |||
|
357 | # Calcultate field components | |||
|
358 | matrix0 = numpy.dot(gdpdtheta,cosmphi) | |||
|
359 | matrix1 = numpy.dot(hdpdtheta,sinmphi) | |||
|
360 | bn[ii] = numpy.dot(rad,(matrix0 + matrix1)) | |||
|
361 | matrix0 = numpy.dot(hp,(mvec*cosmphi)) | |||
|
362 | matrix1 = numpy.dot(gp,(mvec*sinmphi)) | |||
|
363 | be[ii] = numpy.dot((-1*rad),((matrix0 - matrix1)/s)) | |||
|
364 | matrix0 = numpy.dot(gp,cosmphi) | |||
|
365 | matrix1 = numpy.dot(hp,sinmphi) | |||
|
366 | bd[ii] = numpy.dot((-1*nvec*rad),(matrix0 + matrix1)) | |||
|
367 | ||||
|
368 | bmod = numpy.sqrt(bn**2. + be**2. + bd**2.) | |||
|
369 | btheta = numpy.arctan(bd/numpy.sqrt(be**2. + bn**2.))*180/numpy.pi | |||
|
370 | balpha = numpy.arctan(be/bn)*180./numpy.pi | |||
|
371 | ||||
|
372 | #bn : north | |||
|
373 | #be : east | |||
|
374 | #bn : radial | |||
|
375 | #bmod : module | |||
|
376 | ||||
|
377 | ||||
|
378 | return bn, be, bd, bmod, btheta, balpha | |||
|
379 | ||||
|
380 | def str2num(self, datum): | |||
|
381 | try: | |||
|
382 | return int(datum) | |||
|
383 | except: | |||
|
384 | try: | |||
|
385 | return float(datum) | |||
|
386 | except: | |||
|
387 | return datum | |||
|
388 | ||||
|
389 | def __readIGRFfile(self, filename): | |||
|
390 | list_years=[] | |||
|
391 | for i in range(1,24): | |||
|
392 | list_years.append(1895.0 + i*5) | |||
|
393 | ||||
|
394 | epochs=list_years | |||
|
395 | epochs.append(epochs[-1]+5) | |||
|
396 | nepochs = numpy.shape(epochs) | |||
|
397 | ||||
|
398 | gg = numpy.zeros((13,14,nepochs[0]),dtype=float) | |||
|
399 | hh = numpy.zeros((13,14,nepochs[0]),dtype=float) | |||
|
400 | ||||
|
401 | coeffs_file=open(filename) | |||
|
402 | lines=coeffs_file.readlines() | |||
|
403 | ||||
|
404 | coeffs_file.close() | |||
|
405 | ||||
|
406 | for line in lines: | |||
|
407 | items = line.split() | |||
|
408 | g_h = items[0] | |||
|
409 | n = self.str2num(items[1]) | |||
|
410 | m = self.str2num(items[2]) | |||
|
411 | ||||
|
412 | coeffs = items[3:] | |||
|
413 | ||||
|
414 | for i in range(len(coeffs)-1): | |||
|
415 | coeffs[i] = self.str2num(coeffs[i]) | |||
|
416 | ||||
|
417 | #coeffs = numpy.array(coeffs) | |||
|
418 | ncoeffs = numpy.shape(coeffs)[0] | |||
|
419 | ||||
|
420 | if g_h == 'g': | |||
|
421 | # print n," g ",m | |||
|
422 | gg[n-1,m,:]=coeffs | |||
|
423 | elif g_h=='h': | |||
|
424 | # print n," h ",m | |||
|
425 | hh[n-1,m,:]=coeffs | |||
|
426 | # else : | |||
|
427 | # continue | |||
|
428 | ||||
|
429 | # Ultimo Reordenamiento para almacenar . | |||
|
430 | gg[:,:,nepochs[0]-1] = gg[:,:,nepochs[0]-2] + 5*gg[:,:,nepochs[0]-1] | |||
|
431 | hh[:,:,nepochs[0]-1] = hh[:,:,nepochs[0]-2] + 5*hh[:,:,nepochs[0]-1] | |||
|
432 | ||||
|
433 | # return numpy.array([gg,hh]) | |||
|
434 | periods = numpy.array(epochs) | |||
|
435 | g = gg | |||
|
436 | h = hh | |||
|
437 | return periods, g, h | |||
|
438 | ||||
|
439 | ||||
|
440 | def __readIGRFcoeff(self,filename="igrf10coeffs.dat"): | |||
|
441 | """ | |||
|
442 | __readIGRFcoeff reads the coefficients from a binary file which is located in the | |||
|
443 | folder "resource." | |||
|
444 | ||||
|
445 | Parameter | |||
|
446 | --------- | |||
|
447 | filename = A string to specify the name of the file which contains thec coeffs. The | |||
|
448 | default value is "igrf10coeffs.dat" | |||
|
449 | ||||
|
450 | Return | |||
|
451 | ------ | |||
|
452 | periods = A lineal array giving... | |||
|
453 | g1 = | |||
|
454 | h1 = | |||
|
455 | ||||
|
456 | Modification History | |||
|
457 | -------------------- | |||
|
458 | Converted to Python by Freddy R. Galindo, ROJ, 03 October 2009. | |||
|
459 | """ | |||
|
460 | ||||
|
461 | # # igrfile = sys.path[-1] + os.sep + "resource" + os.sep + filename | |||
|
462 | # igrfile = os.path.join('./resource',filename) | |||
|
463 | # f = open(igrfile,'rb') | |||
|
464 | # #f = open(os.getcwd() + os.sep + "resource" + os.sep + filename,'rb') | |||
|
465 | # | |||
|
466 | # # Reading SkyNoise Power (lineal scale) | |||
|
467 | # periods = numpy.fromfile(f,numpy.dtype([('var','<f4')]),23) | |||
|
468 | # periods = periods['var'] | |||
|
469 | # | |||
|
470 | # g = numpy.fromfile(f,numpy.dtype([('var','<f8')]),23*14*14) | |||
|
471 | # g = g['var'].reshape((14,14,23)).transpose() | |||
|
472 | # | |||
|
473 | # h = numpy.fromfile(f,numpy.dtype([('var','<f8')]),23*14*14) | |||
|
474 | # h = h['var'].reshape((14,14,23)).transpose() | |||
|
475 | # | |||
|
476 | # f.close() | |||
|
477 | base_path = os.path.dirname(os.path.abspath(__file__)) | |||
|
478 | filename = os.path.join(base_path, "igrf11coeffs.txt") | |||
|
479 | ||||
|
480 | period_v, g_v, h_v = self.__readIGRFfile(filename) | |||
|
481 | g2 = numpy.zeros((14,14,24)) | |||
|
482 | h2 = numpy.zeros((14,14,24)) | |||
|
483 | g2[1:14,:,:] = g_v | |||
|
484 | h2[1:14,:,:] = h_v | |||
|
485 | ||||
|
486 | g = numpy.transpose(g2, (2,0,1)) | |||
|
487 | h = numpy.transpose(h2, (2,0,1)) | |||
|
488 | periods = period_v.copy() | |||
|
489 | ||||
|
490 | return periods, g, h | |||
|
491 | ||||
|
492 | def rotvector(self,vector,axis=1,ang=0): | |||
|
493 | """ | |||
|
494 | rotvector function returns the new vector generated rotating the rectagular coords. | |||
|
495 | ||||
|
496 | Parameters | |||
|
497 | ---------- | |||
|
498 | vector = A lineal 3-elements array (x,y,z). | |||
|
499 | axis = A integer to specify the axis used to rotate the coord systems. The default | |||
|
500 | value is 1. | |||
|
501 | axis = 1 -> Around "x" | |||
|
502 | axis = 2 -> Around "y" | |||
|
503 | axis = 3 -> Around "z" | |||
|
504 | ang = Angle of rotation (in radians). The default value is zero. | |||
|
505 | ||||
|
506 | Return | |||
|
507 | ------ | |||
|
508 | rotvector = A lineal array of 3 elements giving the new coordinates. | |||
|
509 | ||||
|
510 | Modification History | |||
|
511 | -------------------- | |||
|
512 | Converted to Python by Freddy R. Galindo, ROJ, 01 October 2009. | |||
|
513 | """ | |||
|
514 | ||||
|
515 | if axis==1: | |||
|
516 | t = [[1,0,0],[0,numpy.cos(ang),numpy.sin(ang)],[0,-numpy.sin(ang),numpy.cos(ang)]] | |||
|
517 | elif axis==2: | |||
|
518 | t = [[numpy.cos(ang),0,-numpy.sin(ang)],[0,1,0],[numpy.sin(ang),0,numpy.cos(ang)]] | |||
|
519 | elif axis==3: | |||
|
520 | t = [[numpy.cos(ang),numpy.sin(ang),0],[-numpy.sin(ang),numpy.cos(ang),0],[0,0,1]] | |||
|
521 | ||||
|
522 | rotvector = numpy.array(numpy.dot(numpy.array(t),numpy.array(vector))) | |||
|
523 | ||||
|
524 | return rotvector | |||
|
525 | ||||
|
526 | class AntPatternPlot: | |||
|
527 | ||||
|
528 | def __init__(self): | |||
|
529 | """ | |||
|
530 | AntPatternPlot creates an object to call methods to plot the antenna pattern. | |||
|
531 | ||||
|
532 | Modification History | |||
|
533 | -------------------- | |||
|
534 | Created by Freddy Galindo, ROJ, 06 October 2009. | |||
|
535 | """ | |||
|
536 | ||||
|
537 | self.fig = Figure(figsize=(8,8), facecolor='white') | |||
|
538 | self.ax = self.fig.add_subplot(111) | |||
|
539 | ||||
|
540 | def contPattern(self,iplot=0,gpath='',filename='',mesg='',amp=None ,x=None ,y=None ,getCut=None,title='', save=False): | |||
|
541 | """ | |||
|
542 | contPattern plots a contour map of the antenna pattern. | |||
|
543 | ||||
|
544 | Parameters | |||
|
545 | ---------- | |||
|
546 | iplot = A integer to specify if the plot is the first, second, ... The default va- | |||
|
547 | lue is 0. | |||
|
548 | ||||
|
549 | Examples | |||
|
550 | -------- | |||
|
551 | >> Over_Jro.JroPattern(pattern=2).contPattern() | |||
|
552 | ||||
|
553 | Modification history | |||
|
554 | -------------------- | |||
|
555 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
556 | """ | |||
|
557 | ||||
|
558 | if getCut == 1: | |||
|
559 | return | |||
|
560 | ||||
|
561 | xmax = numpy.max(x) | |||
|
562 | xmin = numpy.min(x) | |||
|
563 | ymax = numpy.max(y) | |||
|
564 | ymin = numpy.min(y) | |||
|
565 | ||||
|
566 | levels = numpy.array([1e-3,1e-2,1e-1,0.5,1.0]) | |||
|
567 | tmp = numpy.round(10*numpy.log10(levels),decimals=1) | |||
|
568 | labels = list(range(5)) | |||
|
569 | for i in numpy.arange(5): | |||
|
570 | labels[i] = str(numpy.int(tmp[i])) | |||
|
571 | ||||
|
572 | ||||
|
573 | colors = ((0,0,1.),(0,170/255.,0),(127/255.,1.,0),(1.,109/255.,0),(128/255.,0,0)) | |||
|
574 | CS = self.ax.contour(x,y,amp.transpose(),levels,colors=colors) | |||
|
575 | fmt = {} | |||
|
576 | for l,s in zip(CS.levels,labels): | |||
|
577 | fmt[l] = s | |||
|
578 | ||||
|
579 | self.ax.annotate('Ng',xy=(-0.05,1.04),xytext=(0.01,0.962),xycoords='axes fraction',arrowprops=dict(facecolor='black', width=1.,shrink=0.2),fontsize=15.) | |||
|
580 | self.ax.annotate(mesg,xy=(0,0),xytext=(0.01,0.01),xycoords='figure fraction') | |||
|
581 | self.ax.clabel(CS,CS.levels,inline=True,fmt=fmt,fontsize=10) | |||
|
582 | self.ax.set_xlim(xmin,xmax) | |||
|
583 | self.ax.set_ylim(ymin,ymax) | |||
|
584 | self.ax.set_title("Total Pattern: " + title) | |||
|
585 | self.ax.set_xlabel("West to South") | |||
|
586 | self.ax.set_ylabel("West to North") | |||
|
587 | self.ax.grid(True) | |||
|
588 | ||||
|
589 | ||||
|
590 | def plotRaDec(self,gpath=None,filename=None,jd=2452640.5,ra_obs=None,xg=None,yg=None,x=None,y=None, save=True): | |||
|
591 | """ | |||
|
592 | plotRaDec draws right ascension and declination lines on a JRO plane. This function | |||
|
593 | must call after conPattern. | |||
|
594 | ||||
|
595 | Parameters | |||
|
596 | ---------- | |||
|
597 | jd = A scalar giving the Julian date. | |||
|
598 | ra_obs = Scalar giving the right ascension of the observatory. | |||
|
599 | xg = A 3-element array to specify .. | |||
|
600 | yg = A 3-element array to specify .. | |||
|
601 | ||||
|
602 | Examples | |||
|
603 | -------- | |||
|
604 | >> Over_Jro.JroPattern(pattern=2).contPattern() | |||
|
605 | >> Over_Jro.JroPattern(pattern=2).plotRaDec(jd=jd,ra_obs=ra_obs,xg=xg,yg=yg) | |||
|
606 | ||||
|
607 | Modification history | |||
|
608 | -------------------- | |||
|
609 | Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009. | |||
|
610 | """ | |||
|
611 | ||||
|
612 | # Finding RA of observatory for a specific date | |||
|
613 | if ra_obs is None:ra_obs = numpy.array([23.37060849]) | |||
|
614 | if xg is None:xg = numpy.array([0.62918474,-0.77725579,0.]) | |||
|
615 | if yg is None:yg = numpy.array([0.77700346,0.62898048,0.02547905]) | |||
|
616 | ||||
|
617 | # Getting HA and DEC axes | |||
|
618 | mindec = -28; maxdec = 4; incdec = 2. | |||
|
619 | ndec = numpy.int((maxdec - mindec)/incdec) + 1 | |||
|
620 | ||||
|
621 | minha = -20; maxha = 20; incha = 2. | |||
|
622 | nha = numpy.int((maxha - minha)/incha) + 1 | |||
|
623 | ||||
|
624 | #mcosx = numpy.zeros((nha,ndec)) | |||
|
625 | #mcosy = numpy.zeros((nha,ndec)) | |||
|
626 | ||||
|
627 | ha_axes = numpy.reshape(numpy.arange(nha)*incha + minha,(nha,1)) | |||
|
628 | ones_dec = numpy.reshape(numpy.zeros(ndec) + 1,(ndec,1)) | |||
|
629 | ha_axes = numpy.dot(ha_axes,ones_dec.transpose()) | |||
|
630 | ha_axes2 = numpy.array(ra_obs - ha_axes) | |||
|
631 | ||||
|
632 | dec_axes = numpy.reshape(numpy.arange(ndec)*incdec + mindec,(ndec,1)) | |||
|
633 | ones_ra = numpy.reshape(numpy.zeros(nha) + 1,(nha,1)) | |||
|
634 | dec_axes = numpy.dot(ones_ra,dec_axes.transpose()) | |||
|
635 | dec_axes2 = numpy.array(dec_axes) | |||
|
636 | ||||
|
637 | ObjHor = Astro_Coords.Equatorial(ha_axes2,dec_axes2,jd) | |||
|
638 | [alt,az,ha] = ObjHor.change2AltAz() | |||
|
639 | ||||
|
640 | z = numpy.transpose(alt)*Misc_Routines.CoFactors.d2r ; z = z.flatten() | |||
|
641 | az = numpy.transpose(az)*Misc_Routines.CoFactors.d2r ; az = az.flatten() | |||
|
642 | ||||
|
643 | vect = numpy.array([numpy.cos(z)*numpy.sin(az),numpy.cos(z)*numpy.cos(az),numpy.sin(z)]) | |||
|
644 | ||||
|
645 | xg = numpy.atleast_2d(xg) | |||
|
646 | dcosx = numpy.array(numpy.dot(xg,vect)) | |||
|
647 | yg = numpy.atleast_2d(yg) | |||
|
648 | dcosy = numpy.array(numpy.dot(yg,vect)) | |||
|
649 | ||||
|
650 | mcosx = dcosx.reshape(ndec,nha) | |||
|
651 | mcosy = dcosy.reshape(ndec,nha) | |||
|
652 | ||||
|
653 | # Defining NAN for points outof limits. | |||
|
654 | xmax = numpy.max(x) | |||
|
655 | xmin = numpy.min(x) | |||
|
656 | ymax = numpy.max(y) | |||
|
657 | ymin = numpy.min(y) | |||
|
658 | ||||
|
659 | factor = 1.3 | |||
|
660 | noval = numpy.where((mcosx>(xmax*factor)) | (mcosx<(xmin*factor))) | |||
|
661 | if noval[0].size>0:mcosx[noval] = numpy.nan | |||
|
662 | noval = numpy.where((mcosy>(ymax*factor)) | (mcosy<(ymin*factor))) | |||
|
663 | if noval[0].size>0:mcosy[noval] = numpy.nan | |||
|
664 | ||||
|
665 | # Plotting HA and declination grid. | |||
|
666 | iha0 = numpy.int((0 - minha)/incha) | |||
|
667 | idec0 = numpy.int((-14 - mindec)/incdec) | |||
|
668 | ||||
|
669 | colorgrid = (1.,109/255.,0) | |||
|
670 | self.ax.plot(mcosx.transpose(),mcosy.transpose(),color=colorgrid,linestyle='--', lw=0.5) | |||
|
671 | for idec in numpy.arange(ndec): | |||
|
672 | if idec != idec0: | |||
|
673 | valx = (mcosx[idec,iha0]<=xmax) & (mcosx[idec,iha0]>=xmin) | |||
|
674 | valy = (mcosy[idec,iha0]<=ymax) & (mcosy[idec,iha0]>=ymin) | |||
|
675 | if valx & valy: | |||
|
676 | text = str(numpy.int(mindec + incdec*idec))+'$^o$' | |||
|
677 | self.ax.text(mcosx[idec,iha0],mcosy[idec,iha0],text) | |||
|
678 | ||||
|
679 | self.ax.plot(mcosx,mcosy,color=colorgrid,linestyle='--',lw=0.5) | |||
|
680 | for iha in numpy.arange(nha): | |||
|
681 | if iha != iha0: | |||
|
682 | valx = (mcosx[idec0,iha]<=xmax) & (mcosx[idec0,iha]>=xmin) | |||
|
683 | valy = (mcosy[idec0,iha]<=ymax) & (mcosy[idec0,iha]>=ymin) | |||
|
684 | if valx & valy: | |||
|
685 | text = str(4*numpy.int(minha + incha*iha))+"'" | |||
|
686 | self.ax.text(mcosx[idec0,iha],mcosy[idec0,iha],text) | |||
|
687 | ||||
|
688 | if save: | |||
|
689 | save_fig = os.path.join(gpath,filename) | |||
|
690 | self.fig.savefig(save_fig,format='png') | |||
|
691 | ||||
|
692 | ||||
|
693 | def plotBField(self,gpath,filename,dcos,alpha, nlon, nlat, dcosxrange, dcosyrange, heights, alpha_i, save=False): | |||
|
694 | """ | |||
|
695 | plotBField draws the magnetic field in a directional cosines plot. | |||
|
696 | ||||
|
697 | Parameters | |||
|
698 | ---------- | |||
|
699 | dcos = An 4-dimensional array giving the directional cosines of the magnetic field | |||
|
700 | over the desired place. | |||
|
701 | alpha = An 3-dimensional array giving the angle of the magnetic field over the desi- | |||
|
702 | red place. | |||
|
703 | nlon = An integer to specify the number of elements per longitude. | |||
|
704 | nlat = An integer to specify the number of elements per latitude. | |||
|
705 | dcosxrange = A 2-element array giving the range of the directional cosines in the | |||
|
706 | "x" axis. | |||
|
707 | dcosyrange = A 2-element array giving the range of the directional cosines in the | |||
|
708 | "y" axis. | |||
|
709 | heights = An array giving the heights (km) where the magnetic field will be modeled By default the magnetic field will be computed at 100, 500 and 1000km. | |||
|
710 | alpha_i = Angle to interpolate the magnetic field. | |||
|
711 | Modification History | |||
|
712 | -------------------- | |||
|
713 | Converted to Python by Freddy R. Galindo, ROJ, 07 October 2009. | |||
|
714 | """ | |||
|
715 | ||||
|
716 | handles = [] | |||
|
717 | objects = [] | |||
|
718 | colors = ['k','m','c','b','g','r','y'] | |||
|
719 | marker = ['-+','-*','-D','-x','-s','->','-o','-^'] | |||
|
720 | ||||
|
721 | alpha_location = numpy.zeros((nlon,2,heights.size)) | |||
|
722 | ||||
|
723 | for ih in numpy.arange(heights.size): | |||
|
724 | alpha_location[:,0,ih] = dcos[:,0,ih,0] | |||
|
725 | for ilon in numpy.arange(nlon): | |||
|
726 | myx = (alpha[ilon,:,ih])[::-1] | |||
|
727 | myy = (dcos[ilon,:,ih,0])[::-1] | |||
|
728 | tck = scipy.interpolate.splrep(myx,myy,s=0) | |||
|
729 | mydcosx = scipy.interpolate.splev(alpha_i,tck,der=0) | |||
|
730 | ||||
|
731 | myx = (alpha[ilon,:,ih])[::-1] | |||
|
732 | myy = (dcos[ilon,:,ih,1])[::-1] | |||
|
733 | tck = scipy.interpolate.splrep(myx,myy,s=0) | |||
|
734 | mydcosy = scipy.interpolate.splev(alpha_i,tck,der=0) | |||
|
735 | alpha_location[ilon,:,ih] = numpy.array([mydcosx, mydcosy]) | |||
|
736 | ||||
|
737 | ||||
|
738 | ObjFig, = self.ax.plot(alpha_location[:,0,ih],alpha_location[:,1,ih], | |||
|
739 | marker[ih % 8],color=colors[numpy.int(ih/8)],ms=4.5,lw=0.5) | |||
|
740 | handles.append(ObjFig) | |||
|
741 | objects.append(numpy.str(heights[ih]) + ' km') | |||
|
742 | ||||
|
743 | legend = self.ax.legend(handles, objects,loc="lower right", numpoints=1, handlelength=0.3, | |||
|
744 | handletextpad=0.02, borderpad=0.3, labelspacing=0.1) | |||
|
745 | self.ax.add_artist(legend) | |||
|
746 | ||||
|
747 | def plotCelestial(self, jd, main_dec, tod, maxha_min, objects, glat, glon, xg, yg, dcosxrange, dcosyrange): | |||
|
748 | ||||
|
749 | self.tod = tod | |||
|
750 | ||||
|
751 | self.dcosx_sun = 1 | |||
|
752 | self.dcosy_sun = 1 | |||
|
753 | self.ha_sun = 1 | |||
|
754 | self.time_sun = 1 | |||
|
755 | ||||
|
756 | self.dcosx_moon = 1 | |||
|
757 | self.dcosy_moon = 1 | |||
|
758 | self.ha_moon = 1 | |||
|
759 | self.time_moon = 1 | |||
|
760 | ||||
|
761 | self.dcosx_hydra = 1 | |||
|
762 | self.dcosy_hydra = 1 | |||
|
763 | self.ha_hydra = 1 | |||
|
764 | self.time_hydra = 1 | |||
|
765 | ||||
|
766 | self.dcosx_galaxy = 1 | |||
|
767 | self.dcosy_galaxy = 1 | |||
|
768 | self.ha_galaxy = 1 | |||
|
769 | self.time_galaxy = 1 | |||
|
770 | ||||
|
771 | tod = self.tod | |||
|
772 | ||||
|
773 | maxlev = 24; minlev = 0; maxcol = 39; mincol = 10 | |||
|
774 | handles = [] | |||
|
775 | titles = ['$Sun$','$Moon$','$Hydra$','$Galaxy$'] | |||
|
776 | marker = ['--^','--s','--*','--o'] | |||
|
777 | ||||
|
778 | # Getting RGB table to plot celestial object over Jicamarca | |||
|
779 | colortable = ['olive', 'indigo', 'cyan', 'red'] | |||
|
780 | labels = [] | |||
|
781 | for io in objects: | |||
|
782 | ObjBodies = Astro_Coords.CelestialBodies() | |||
|
783 | if io==0: | |||
|
784 | [ra,dec,sunlon,sunobliq] = ObjBodies.sunpos(jd) | |||
|
785 | elif io==1: | |||
|
786 | [ra,dec,dist,moonlon,moonlat] = ObjBodies.moonpos(jd) | |||
|
787 | elif io==2: | |||
|
788 | [ra,dec] = ObjBodies.hydrapos() | |||
|
789 | elif io==3: | |||
|
790 | [maxra,ra] = ObjBodies.skynoise_jro(dec_cut=main_dec) | |||
|
791 | ra = maxra*15. | |||
|
792 | dec = main_dec | |||
|
793 | ||||
|
794 | ObjEq = Astro_Coords.Equatorial(ra, dec, jd, lat=glat, lon=glon) | |||
|
795 | [alt, az, ha] = ObjEq.change2AltAz() | |||
|
796 | vect = numpy.array([az,alt]).transpose() | |||
|
797 | vect = Misc_Routines.Vector(vect,direction=0).Polar2Rect() | |||
|
798 | ||||
|
799 | dcosx = numpy.array(numpy.dot(vect,xg)) | |||
|
800 | dcosy = numpy.array(numpy.dot(vect,yg)) | |||
|
801 | wrap = numpy.where(ha>=180.) | |||
|
802 | ||||
|
803 | if wrap[0].size>0: | |||
|
804 | ha[wrap] = ha[wrap] - 360. | |||
|
805 | ||||
|
806 | val = numpy.where((numpy.abs(ha))<=(maxha_min*0.25)) | |||
|
807 | ||||
|
808 | if val[0].size>2: | |||
|
809 | tod_1 = tod*1. | |||
|
810 | shift_1 = numpy.where(tod>12.) | |||
|
811 | tod_1[shift_1] = tod_1[shift_1] - 24. | |||
|
812 | tod_2 = tod*1. | |||
|
813 | shift_2 = numpy.where(tod<12.) | |||
|
814 | tod_2[shift_2] = tod_2[shift_2] + 24. | |||
|
815 | ||||
|
816 | diff0 = numpy.nanmax(tod[val]) - numpy.nanmin(tod[val]) | |||
|
817 | diff1 = numpy.nanmax(tod_1[val]) - numpy.nanmin(tod_1[val]) | |||
|
818 | diff2 = numpy.nanmax(tod_2[val]) - numpy.nanmin(tod_2[val]) | |||
|
819 | ||||
|
820 | if ((diff0<=diff1) & (diff0<=diff2)): | |||
|
821 | tod_0 = tod | |||
|
822 | elif ((diff1<diff0) & (diff1<diff2)): | |||
|
823 | tod_0 = tod_1 | |||
|
824 | else: | |||
|
825 | tod_0 = tod_2 | |||
|
826 | ||||
|
827 | if io==0: | |||
|
828 | self.dcosx_sun = dcosx[val] | |||
|
829 | self.dcosy_sun = dcosy[val] | |||
|
830 | self.ha_sun = ha[val] | |||
|
831 | self.time_sun = numpy.median(tod_0[val]) | |||
|
832 | elif io==1: | |||
|
833 | self.dcosx_moon = dcosx[val] | |||
|
834 | self.dcosy_moon = dcosy[val] | |||
|
835 | self.ha_moon = ha[val] | |||
|
836 | self.time_moon = numpy.median(tod_0[val]) | |||
|
837 | elif io==2: | |||
|
838 | self.dcosx_hydra = dcosx[val] | |||
|
839 | self.dcosy_hydra = dcosy[val] | |||
|
840 | self.ha_hydra = ha[val] | |||
|
841 | self.time_hydra = numpy.mean(tod_0[val]) | |||
|
842 | elif io==3: | |||
|
843 | self.dcosx_galaxy = dcosx[val] | |||
|
844 | self.dcosy_galaxy = dcosy[val] | |||
|
845 | self.ha_galaxy = ha[val] | |||
|
846 | self.time_galaxy = numpy.mean(tod_0[val]) | |||
|
847 | ||||
|
848 | index = numpy.mean(tod_0[val]) - minlev | |||
|
849 | index = (index*(maxcol - mincol)/(maxlev - minlev)) + mincol | |||
|
850 | index = numpy.int(index) | |||
|
851 | figobjects, = self.ax.plot(dcosx[val],dcosy[val],marker[io],\ | |||
|
852 | lw=1,ms=7,mew=0,color=colortable[io]) | |||
|
853 | handles.append(figobjects) | |||
|
854 | labels.append(titles[io]) | |||
|
855 | ||||
|
856 | ||||
|
857 | legend = self.ax.legend(handles,labels,loc="lower left", numpoints=1, handlelength=0.5, \ | |||
|
858 | borderpad=0.3, handletextpad=0.1,labelspacing=0.2,fontsize='small') | |||
|
859 | ||||
|
860 | self.ax.add_artist(legend) | |||
|
861 | ||||
|
862 | ||||
|
863 | class JroPattern(): | |||
|
864 | def __init__(self,pattern=0,path=None,filename=None,nptsx=101,nptsy=101,maxphi=5,fftopt=0, \ | |||
|
865 | getcut=0,dcosx=None,dcosy=None,eomwl=6,airwl=4, **kwargs): | |||
|
866 | """ | |||
|
867 | JroPattern class creates an object to represent the useful parameters for beam mode- | |||
|
868 | lling of the Jicamarca VHF radar. | |||
|
869 | ||||
|
870 | Parameters | |||
|
871 | ---------- | |||
|
872 | pattern = An integer (See JroAntSetup to know the available values) to load a prede- | |||
|
873 | fined configuration. The default value is 0. To use a user-defined configuration | |||
|
874 | pattern must be None. | |||
|
875 | path = A string giving the directory that contains the user-configuration file. PATH | |||
|
876 | will work if pattern is None. | |||
|
877 | filename = A string giving the name of the user-configuration file. FILENAME will | |||
|
878 | work if pattern is None. | |||
|
879 | nptsx = A scalar to specify the number of points used to define the angular resolu- | |||
|
880 | tion in the "x" axis. The default value is 101. | |||
|
881 | nptsy = A scalar to specify the number of points used to define the angular resolu- | |||
|
882 | tion in the "x" axis. The default value is 101. | |||
|
883 | maxphi = A scalar giving the maximum (absolute) angle (in degree) to model the ante- | |||
|
884 | nna pattern. The default value is 5 degrees. | |||
|
885 | fftopt = Set this input to 1 to model the beam using FFT. To model using antenna | |||
|
886 | theory set to 0 (default value). | |||
|
887 | getcut = Set to 1 to show an antenna cut instead of a contour plot of itself (set to | |||
|
888 | 0). The defautl value is 0. | |||
|
889 | dcosx = An array giving the directional cosines for the x-axis. DCOSX will work if | |||
|
890 | getcut is actived. | |||
|
891 | dcosy = An array giving the directional cosines for the y-axis. DCOSY will work if | |||
|
892 | getcut is actived. | |||
|
893 | eomwl = A scalar giving the radar wavelength. The default value is 6m (50 MHZ). | |||
|
894 | airwl = Set this input to float (or intger) to specify the wavelength (in meters) of | |||
|
895 | the transmitted EOM wave in the air. The default value is 4m. | |||
|
896 | ||||
|
897 | Modification History | |||
|
898 | -------------------- | |||
|
899 | Converted to Object-oriented Programming by Freddy Galindo, ROJ, 20 September 2009. | |||
|
900 | """ | |||
|
901 | ||||
|
902 | ||||
|
903 | ||||
|
904 | # Getting antenna configuration. | |||
|
905 | if filename: | |||
|
906 | setup = JroAntSetup.ReturnSetup(path=path,filename=filename,pattern=pattern) | |||
|
907 | ||||
|
908 | ues = setup["ues"] | |||
|
909 | phase = setup["phase"] | |||
|
910 | gaintx = setup["gaintx"] | |||
|
911 | gainrx = setup["gainrx"] | |||
|
912 | justrx = setup["justrx"] | |||
|
913 | self.title = setup["title"] | |||
|
914 | else: | |||
|
915 | ues = kwargs["ues"] | |||
|
916 | phase = kwargs["phases"] | |||
|
917 | gaintx = kwargs["gain_tx"] | |||
|
918 | gainrx = kwargs["gain_rx"] | |||
|
919 | justrx = kwargs["just_rx"] | |||
|
920 | self.title = kwargs.get("title", "JRO Pattern") | |||
|
921 | ||||
|
922 | # Defining attributes for JroPattern class. | |||
|
923 | # Antenna configuration | |||
|
924 | ||||
|
925 | self.uestx = ues | |||
|
926 | self.phasetx = phase | |||
|
927 | self.gaintx = gaintx | |||
|
928 | self.uesrx = ues | |||
|
929 | self.phaserx = phase | |||
|
930 | self.gainrx = gainrx | |||
|
931 | self.justrx = justrx | |||
|
932 | ||||
|
933 | # Pattern resolution & method to model | |||
|
934 | self.maxphi = maxphi | |||
|
935 | self.nptsx = nptsx | |||
|
936 | self.nptsy = nptsy | |||
|
937 | self.fftopt = fftopt | |||
|
938 | ||||
|
939 | # To get a cut of the pattern. | |||
|
940 | self.getcut = getcut | |||
|
941 | ||||
|
942 | maxdcos = numpy.sin(maxphi*Misc_Routines.CoFactors.d2r) | |||
|
943 | if dcosx==None:dcosx = ((numpy.arange(nptsx,dtype=float)/(nptsx-1))-0.5)*2*maxdcos | |||
|
944 | if dcosy==None:dcosy = ((numpy.arange(nptsy,dtype=float)/(nptsy-1))-0.5)*2*maxdcos | |||
|
945 | self.dcosx = dcosx | |||
|
946 | self.dcosy = dcosy | |||
|
947 | self.nx = dcosx.size | |||
|
948 | self.ny = dcosy.size*(getcut==0) + (getcut==1) | |||
|
949 | ||||
|
950 | self.eomwl = eomwl | |||
|
951 | self.airwl = airwl | |||
|
952 | ||||
|
953 | self.kk = 2.*numpy.pi/eomwl | |||
|
954 | ||||
|
955 | self.pattern = None | |||
|
956 | self.meanpos = None | |||
|
957 | self.norpattern = None | |||
|
958 | self.maxpattern = None | |||
|
959 | ||||
|
960 | ||||
|
961 | ||||
|
962 | self.getPattern() | |||
|
963 | ||||
|
964 | def getPattern(self): | |||
|
965 | """ | |||
|
966 | getpattern method returns the modeled total antenna pattern and its mean position. | |||
|
967 | ||||
|
968 | Return | |||
|
969 | ------ | |||
|
970 | pattern = An array giving the Modelled antenna pattern. | |||
|
971 | mean_pos = A 2-elements array giving the mean position of the main beam. | |||
|
972 | ||||
|
973 | Examples | |||
|
974 | -------- | |||
|
975 | >> [pattern, mean_pos] = JroPattern(pattern=2).getPattern() | |||
|
976 | >> print meanpos | |||
|
977 | [ 8.08728085e-14 -4.78193873e-14] | |||
|
978 | ||||
|
979 | Modification history | |||
|
980 | -------------------- | |||
|
981 | Developed by Jorge L. Chau. | |||
|
982 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
983 | """ | |||
|
984 | ||||
|
985 | if (self.fftopt>0) and (self.getcut>0): | |||
|
986 | #print "Conflict bewteen fftopt and getcut" | |||
|
987 | #print "To get a cut of the antenna pattern uses ffopt=0" | |||
|
988 | return None, None | |||
|
989 | ||||
|
990 | if (self.fftopt==0): | |||
|
991 | # Getting antenna pattern using the array method | |||
|
992 | self.pattern = self.__usingArray(rx=1) | |||
|
993 | if (self.justrx==0):self.pattern = self.pattern*self.__usingArray(rx=0) | |||
|
994 | ||||
|
995 | elif (self.fftopt>0): | |||
|
996 | # Getting antenna pattern using FFT method | |||
|
997 | self.pattern = self.__usingFFT(rx=1) | |||
|
998 | if (self.justrx==0):self.pattern = self.pattern*self.__usingFFT(rx=0) | |||
|
999 | ||||
|
1000 | self.maxpattern = numpy.nanmax(self.pattern) | |||
|
1001 | self.norpattern = self.pattern/self.maxpattern | |||
|
1002 | if self.getcut==0:self.__getBeamPars() | |||
|
1003 | ||||
|
1004 | def __usingArray(self,rx): | |||
|
1005 | """ | |||
|
1006 | __usingArray method returns the Jicamarca antenna pattern computed using array model | |||
|
1007 | ||||
|
1008 | pattern = dipolepattern x modulepattern | |||
|
1009 | ||||
|
1010 | Parameters | |||
|
1011 | ---------- | |||
|
1012 | rx = Set to 1 to use the Rx information. Otherwise set to 0 for Tx. | |||
|
1013 | ||||
|
1014 | Return | |||
|
1015 | ------ | |||
|
1016 | pattern = An array giving the modelled antenna pattern using the array model. | |||
|
1017 | ||||
|
1018 | Modification history | |||
|
1019 | -------------------- | |||
|
1020 | Developed by Jorge L. Chau. | |||
|
1021 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
1022 | """ | |||
|
1023 | ||||
|
1024 | if rx==1: | |||
|
1025 | ues = self.uesrx | |||
|
1026 | phase = self.phaserx | |||
|
1027 | gain = self.gainrx | |||
|
1028 | elif rx==0: | |||
|
1029 | ues = self.uestx | |||
|
1030 | phase = self.phasetx | |||
|
1031 | gain = self.gaintx | |||
|
1032 | ||||
|
1033 | ues = ues*360./self.airwl | |||
|
1034 | phase = phase*360./self.airwl | |||
|
1035 | ||||
|
1036 | for ii in range(4): | |||
|
1037 | if ii==0:dim = numpy.array([4,0,8,4]) # WEST | |||
|
1038 | elif ii==1:dim = numpy.array([0,0,4,4]) # NORTH | |||
|
1039 | elif ii==2:dim = numpy.array([0,4,4,8]) # EAST | |||
|
1040 | elif ii==3:dim = numpy.array([4,4,8,8]) # SOUTH | |||
|
1041 | xi = dim[0]; xf = dim[2]; yi = dim[1]; yf = dim[3] | |||
|
1042 | phase[xi:xf,yi:yf] = phase[xi:xf,yi:yf] + ues[ii] | |||
|
1043 | ||||
|
1044 | phase = -phase | |||
|
1045 | ||||
|
1046 | ar = self.eomwl*numpy.array([[0.5,6., 24.5],[0.5,6.,24.5]]) | |||
|
1047 | nr = numpy.array([[12.,4.,2.],[12.,4.,2.]]) | |||
|
1048 | lr = 0.25*self.eomwl*numpy.array([[0,0.,0],[0.,0,0]]) | |||
|
1049 | ||||
|
1050 | # Computing module and dipole patterns. | |||
|
1051 | pattern = (numpy.abs(self.__dipPattern(ar,nr,lr)*self.__modPattern(phase,gain)))**2 | |||
|
1052 | ||||
|
1053 | return pattern | |||
|
1054 | ||||
|
1055 | def __usingFFT(self,rx): | |||
|
1056 | """ | |||
|
1057 | __usingFFT method returns the Jicamarca antenna pattern computed using The Fast Fou- | |||
|
1058 | rier Transform. | |||
|
1059 | ||||
|
1060 | pattern = iFFT(FFT(gain*EXP(j*phase))) | |||
|
1061 | ||||
|
1062 | Parameters | |||
|
1063 | ---------- | |||
|
1064 | rx = Set to 1 to use the Rx information. Otherwise set to 0 for Tx. | |||
|
1065 | ||||
|
1066 | Return | |||
|
1067 | ------ | |||
|
1068 | pattern = An array giving the modelled antenna pattern using the array model. | |||
|
1069 | ||||
|
1070 | Modification history | |||
|
1071 | -------------------- | |||
|
1072 | Developed by Jorge L. Chau. | |||
|
1073 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
1074 | """ | |||
|
1075 | ||||
|
1076 | if rx==1: | |||
|
1077 | ues = self.uesrx | |||
|
1078 | phase = self.phaserx | |||
|
1079 | gain = self.gainrx | |||
|
1080 | elif rx==0: | |||
|
1081 | ues = self.uestx | |||
|
1082 | phase = self.phasetx | |||
|
1083 | gain = self.gaintx | |||
|
1084 | ||||
|
1085 | ues = ues*360./self.airwl | |||
|
1086 | phase = phase*360./self.airwl | |||
|
1087 | ||||
|
1088 | for ii in range(4): | |||
|
1089 | if ii==0:dim = numpy.array([4,0,8,4]) # WEST | |||
|
1090 | elif ii==1:dim = numpy.array([0,0,4,4]) # NORTH | |||
|
1091 | elif ii==2:dim = numpy.array([0,4,4,8]) # EAST | |||
|
1092 | elif ii==3:dim = numpy.array([4,4,8,8]) # SOUTH | |||
|
1093 | xi = dim[0]; xf = dim[2]; yi = dim[1]; yf = dim[3] | |||
|
1094 | phase[xi:xf,yi:yf] = phase[xi:xf,yi:yf] + ues[ii] | |||
|
1095 | ||||
|
1096 | phase = -phase | |||
|
1097 | ||||
|
1098 | delta_x = self.eomwl/2. | |||
|
1099 | delta_y = self.eomwl/2. | |||
|
1100 | ||||
|
1101 | nxfft = 2048 | |||
|
1102 | nyfft = 2048 | |||
|
1103 | dcosx = (numpy.arange(nxfft) - (0.5*nxfft))/(nxfft*delta_x)*self.eomwl | |||
|
1104 | dcosy = (numpy.arange(nyfft) - (0.5*nyfft))/(nyfft*delta_y)*self.eomwl | |||
|
1105 | ||||
|
1106 | fft_gain = numpy.zeros((nxfft,nyfft)) | |||
|
1107 | fft_phase = numpy.zeros((nxfft,nyfft)) | |||
|
1108 | ||||
|
1109 | nx = 8 | |||
|
1110 | ny = 8 | |||
|
1111 | ndx =12 | |||
|
1112 | ndy =12 | |||
|
1113 | for iy in numpy.arange(ny): | |||
|
1114 | for ix in numpy.arange(nx): | |||
|
1115 | ix1 = nxfft/2-self.nx/2*ndx+ix*ndx | |||
|
1116 | if ix<(nx/2):ix1 = ix1 - 1 | |||
|
1117 | if ix>=(nx/2):ix1 = ix1 + 1 | |||
|
1118 | ||||
|
1119 | iy1 = nyfft/2-ny/2*ndx+iy*ndy | |||
|
1120 | if iy<(ny/2):iy1 = iy1 - 1 | |||
|
1121 | if iy>=(ny/2):iy1 = iy1 + 1 | |||
|
1122 | ||||
|
1123 | fft_gain[ix1:ix1+ndx-1,iy1:iy1+ndy-1] = gain[ix,ny-1-iy] | |||
|
1124 | fft_phase[ix1:ix1+ndx-1,iy1:iy1+ndy-1] = phase[ix,ny-1-iy] | |||
|
1125 | ||||
|
1126 | ||||
|
1127 | fft_phase = fft_phase*Misc_Routines.CoFactors.d2r | |||
|
1128 | ||||
|
1129 | pattern = numpy.abs(numpy.fft.fft2(fft_gain*numpy.exp(numpy.complex(0,1)*fft_phase)))**2 | |||
|
1130 | pattern = numpy.fft.fftshift(pattern) | |||
|
1131 | ||||
|
1132 | xvals = numpy.where((dcosx>=(numpy.min(self.dcosx))) & (dcosx<=(numpy.max(self.dcosx)))) | |||
|
1133 | yvals = numpy.where((dcosy>=(numpy.min(self.dcosy))) & (dcosy<=(numpy.max(self.dcosy)))) | |||
|
1134 | ||||
|
1135 | pattern = pattern[xvals[0][0]:xvals[0][-1],yvals[0][0]:yvals[0][-1]] | |||
|
1136 | ||||
|
1137 | return pattern | |||
|
1138 | ||||
|
1139 | ||||
|
1140 | def __dipPattern(self,ar,nr,lr): | |||
|
1141 | """ | |||
|
1142 | _dipPattern function computes the dipole's pattern to the Jicamarca radar. The next | |||
|
1143 | equation defines the pattern as a function of the mainlobe direction: | |||
|
1144 | ||||
|
1145 | sincx = SIN(k/2*n0x*(a0x*SIN(phi)*COS(alpha)))/SIN(k/2*(a0x*SIN(phi)*COS(alpha))) | |||
|
1146 | sincy = SIN(k/2*n0y*(a0y*SIN(phi)*SIN(alpha)))/SIN(k/2*(a0y*SIN(phi)*SIN(alpha))) | |||
|
1147 | A0(phi,alpha) = sincx*sincy | |||
|
1148 | Parameters | |||
|
1149 | ---------- | |||
|
1150 | ar = ? | |||
|
1151 | nr = ? | |||
|
1152 | lr = ? | |||
|
1153 | ||||
|
1154 | Return | |||
|
1155 | ------ | |||
|
1156 | dipole = An array giving antenna pattern from the dipole point of view.. | |||
|
1157 | ||||
|
1158 | Modification history | |||
|
1159 | -------------------- | |||
|
1160 | Developed by Jorge L. Chau. | |||
|
1161 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
1162 | """ | |||
|
1163 | ||||
|
1164 | dipole = numpy.zeros((self.nx,self.ny),dtype=complex) | |||
|
1165 | for iy in range(self.ny): | |||
|
1166 | for ix in range(self.nx): | |||
|
1167 | yindex = iy*(self.getcut==0) + ix*(self.getcut==1) | |||
|
1168 | ||||
|
1169 | argx = ar[0,0]*self.dcosx[ix] - lr[0,0] | |||
|
1170 | if argx == 0.0: | |||
|
1171 | junkx = nr[0,0] | |||
|
1172 | else: | |||
|
1173 | junkx = numpy.sin(0.5*self.kk*nr[0,0]*argx)/numpy.sin(0.5*self.kk*argx) | |||
|
1174 | ||||
|
1175 | ||||
|
1176 | argy = ar[1,0]*self.dcosy[yindex] - lr[1,0] | |||
|
1177 | if argy == 0.0: | |||
|
1178 | junky = nr[1,0] | |||
|
1179 | else: | |||
|
1180 | junky = numpy.sin(0.5*self.kk*nr[1,0]*argy)/numpy.sin(0.5*self.kk*argy) | |||
|
1181 | ||||
|
1182 | ||||
|
1183 | dipole[ix,iy] = junkx*junky | |||
|
1184 | ||||
|
1185 | return dipole | |||
|
1186 | ||||
|
1187 | def __modPattern(self,phase,gain): | |||
|
1188 | """ | |||
|
1189 | ModPattern computes the module's pattern to the Jicamarca radar. The next equation | |||
|
1190 | defines the pattern as a function mainlobe direction: | |||
|
1191 | ||||
|
1192 | phasex = pos(x)*SIN(phi)*COS(alpha) | |||
|
1193 | phasey = pos(y)*SIN(phi)*SIN(alpha) | |||
|
1194 | ||||
|
1195 | A1(phi,alpha) = TOTAL(gain*EXP(COMPLEX(0,k*(phasex+phasey)+phase))) | |||
|
1196 | ||||
|
1197 | Parameters | |||
|
1198 | ---------- | |||
|
1199 | phase = Bidimensional array (8x8) giving the phase (in meters) of each module. | |||
|
1200 | gain = Bidimensional array (8x8) giving to define modules will be active (ones) | |||
|
1201 | and which will not (zeros). | |||
|
1202 | ||||
|
1203 | Return | |||
|
1204 | ------ | |||
|
1205 | module = An array giving antenna pattern from the module point of view.. | |||
|
1206 | ||||
|
1207 | Modification history | |||
|
1208 | -------------------- | |||
|
1209 | Developed by Jorge L. Chau. | |||
|
1210 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
1211 | """ | |||
|
1212 | ||||
|
1213 | pos = self.eomwl*attenuation | |||
|
1214 | posx = pos[0,:,:] | |||
|
1215 | posy = pos[1,:,:] | |||
|
1216 | ||||
|
1217 | phase = phase*Misc_Routines.CoFactors.d2r | |||
|
1218 | module = numpy.zeros((self.nx,self.ny),dtype=complex) | |||
|
1219 | for iy in range(self.ny): | |||
|
1220 | for ix in range(self.nx): | |||
|
1221 | yindex = iy*(self.getcut==0) + ix*(self.getcut==1) | |||
|
1222 | phasex = posx*self.dcosx[ix] | |||
|
1223 | phasey = posy*self.dcosy[yindex] | |||
|
1224 | tmp = gain*numpy.exp(numpy.complex(0,1.)*(self.kk*(phasex+phasey)+phase)) | |||
|
1225 | module[ix,iy] = tmp.sum() | |||
|
1226 | ||||
|
1227 | return module | |||
|
1228 | ||||
|
1229 | def __getBeamPars(self): | |||
|
1230 | """ | |||
|
1231 | _getBeamPars computes the main-beam parameters of the antenna. | |||
|
1232 | ||||
|
1233 | Modification history | |||
|
1234 | -------------------- | |||
|
1235 | Developed by Jorge L. Chau. | |||
|
1236 | Converted to Python by Freddy R. Galindo, ROJ, 20 September 2009. | |||
|
1237 | """ | |||
|
1238 | ||||
|
1239 | dx = self.dcosx[1] - self.dcosx[0] | |||
|
1240 | dy = self.dcosy[1] - self.dcosy[0] | |||
|
1241 | ||||
|
1242 | amp = self.norpattern | |||
|
1243 | ||||
|
1244 | xx = numpy.resize(self.dcosx,(self.nx,self.nx)).transpose() | |||
|
1245 | yy = numpy.resize(self.dcosy,(self.ny,self.ny)) | |||
|
1246 | ||||
|
1247 | mm0 = amp[numpy.where(amp > 0.5)] | |||
|
1248 | xx0 = xx[numpy.where(amp > 0.5)] | |||
|
1249 | yy0 = yy[numpy.where(amp > 0.5)] | |||
|
1250 | ||||
|
1251 | xc = numpy.sum(mm0*xx0)/numpy.sum(mm0) | |||
|
1252 | yc = numpy.sum(mm0*yy0)/numpy.sum(mm0) | |||
|
1253 | rc = numpy.sqrt(mm0.size*dx*dy/numpy.pi) | |||
|
1254 | ||||
|
1255 | nnx = numpy.where(numpy.abs(self.dcosx - xc) < rc) | |||
|
1256 | nny = numpy.where(numpy.abs(self.dcosy - yc) < rc) | |||
|
1257 | ||||
|
1258 | mm1 = amp[numpy.min(nnx):numpy.max(nnx)+1,numpy.min(nny):numpy.max(nny)+1] | |||
|
1259 | xx1 = self.dcosx[numpy.min(nnx):numpy.max(nnx)+1] | |||
|
1260 | yy1 = self.dcosy[numpy.min(nny):numpy.max(nny)+1] | |||
|
1261 | ||||
|
1262 | # fitting data into the main beam. | |||
|
1263 | ||||
|
1264 | params = gaussfit.fitgaussian(mm1) | |||
|
1265 | ||||
|
1266 | # Tranforming from indexes to axis' values | |||
|
1267 | xcenter = xx1[0] + (((xx1[xx1.size-1] - xx1[0])/(xx1.size -1))*(params[1])) | |||
|
1268 | ycenter = yy1[0] + (((yy1[yy1.size-1] - yy1[0])/(yy1.size -1))*(params[2])) | |||
|
1269 | xwidth = ((xx1[xx1.size-1] - xx1[0])/(xx1.size-1))*(params[3])*(1/Misc_Routines.CoFactors.d2r) | |||
|
1270 | ywidth = ((yy1[yy1.size-1] - yy1[0])/(yy1.size-1))*(params[4])*(1/Misc_Routines.CoFactors.d2r) | |||
|
1271 | meanwx = (xwidth*ywidth) | |||
|
1272 | meanpos = numpy.array([xcenter,ycenter]) | |||
|
1273 | ||||
|
1274 | #print 'Position: %f %f' %(xcenter,ycenter) | |||
|
1275 | #print 'Widths: %f %f' %(xwidth, ywidth) | |||
|
1276 | #print 'BWHP: %f' %(2*numpy.sqrt(2*meanwx)*numpy.sqrt(-numpy.log(0.5))) | |||
|
1277 | ||||
|
1278 | self.meanpos = meanpos | |||
|
1279 | ||||
|
1280 | ||||
|
1281 | class overJroShow: | |||
|
1282 | ||||
|
1283 | __serverdocspath = '' | |||
|
1284 | __tmpDir = '' | |||
|
1285 | ||||
|
1286 | def __init__(self, title='', heights=None): | |||
|
1287 | self.year = None | |||
|
1288 | self.month = None | |||
|
1289 | self.dom = None | |||
|
1290 | self.pattern = None | |||
|
1291 | self.maxphi = None | |||
|
1292 | self.heights = None | |||
|
1293 | self.filename = None | |||
|
1294 | self.showType = None | |||
|
1295 | self.path = None | |||
|
1296 | self.objects = None | |||
|
1297 | self.nptsx = 101 | |||
|
1298 | self.nptsy = 101 | |||
|
1299 | self.fftopt = 0 | |||
|
1300 | self.site = 1 | |||
|
1301 | self.dcosx = 1 | |||
|
1302 | self.dcosy = 1 | |||
|
1303 | self.dcosxrange = None | |||
|
1304 | self.dcosyrange = None | |||
|
1305 | self.maxha_min= 0. | |||
|
1306 | self.show_object = None | |||
|
1307 | self.dcosx_mag = None | |||
|
1308 | self.dcosy_mag = None | |||
|
1309 | self.ha_mag = None | |||
|
1310 | self.time_mag = None | |||
|
1311 | self.main_dec = None | |||
|
1312 | self.ObjC = None | |||
|
1313 | self.ptitle = title | |||
|
1314 | self.path4plotname = None | |||
|
1315 | self.plotname0 = None | |||
|
1316 | self.plotname1 = None | |||
|
1317 | self.plotname2 = None | |||
|
1318 | self.scriptHeaders = 0 | |||
|
1319 | self.glat = -11.95 | |||
|
1320 | self.glon = -76.8667 | |||
|
1321 | self.UT = 5 #timezone | |||
|
1322 | ||||
|
1323 | self.glat = -11.951481 | |||
|
1324 | self.glon = -76.874383 | |||
|
1325 | if heights is None: | |||
|
1326 | self.heights = numpy.array([100.,500.,1000.]) | |||
|
1327 | else: | |||
|
1328 | self.heights = numpy.array(heights) | |||
|
1329 | ||||
|
1330 | def initParameters1(self): | |||
|
1331 | ||||
|
1332 | gui=1 | |||
|
1333 | if self.pattern==None: | |||
|
1334 | if gui==1: self.filename = self.filename.split(',') | |||
|
1335 | ||||
|
1336 | pattern = numpy.atleast_1d(self.pattern) | |||
|
1337 | filename = numpy.atleast_1d(self.filename) | |||
|
1338 | ||||
|
1339 | npatterns = numpy.max(numpy.array([pattern.size,filename.size])) | |||
|
1340 | ||||
|
1341 | self.pattern = numpy.resize(pattern,npatterns) | |||
|
1342 | self.filename = numpy.resize(filename,npatterns) | |||
|
1343 | ||||
|
1344 | self.doy = datetime.datetime(self.year,self.month,self.dom).timetuple().tm_yday | |||
|
1345 | ||||
|
1346 | ||||
|
1347 | if self.objects==None: | |||
|
1348 | self.objects=numpy.zeros(5) | |||
|
1349 | else: | |||
|
1350 | tmp = numpy.atleast_1d(self.objects) | |||
|
1351 | self.objects = numpy.zeros(5) | |||
|
1352 | self.objects[0:tmp.size] = tmp | |||
|
1353 | ||||
|
1354 | self.show_object = self.objects | |||
|
1355 | ||||
|
1356 | self.maxha_min = 4*self.maxphi*numpy.sqrt(2)*1.25 | |||
|
1357 | ||||
|
1358 | ||||
|
1359 | ||||
|
1360 | #ROJ geographic coordinates and time zone | |||
|
1361 | self.UT = 5 #timezone | |||
|
1362 | self.glat = -11.951481 | |||
|
1363 | self.glon = -76.874383 | |||
|
1364 | ||||
|
1365 | ||||
|
1366 | self.junkjd = TimeTools.Time(self.year,self.month,self.dom).change2julday() | |||
|
1367 | self.junklst = TimeTools.Julian(self.junkjd).change2lst(longitude=self.glon) | |||
|
1368 | ||||
|
1369 | # Finding RA of observatory for a specific date | |||
|
1370 | self.ra_obs = self.junklst*Misc_Routines.CoFactors.h2d | |||
|
1371 | ||||
|
1372 | def initParameters(self, dt): | |||
|
1373 | ||||
|
1374 | self.year = dt.year | |||
|
1375 | self.month = dt.month | |||
|
1376 | self.dom = dt.day | |||
|
1377 | # Defining plot filenames | |||
|
1378 | self.path4plotname = os.path.join(self.__serverdocspath,self.__tmpDir) | |||
|
1379 | self.plotname0 = 'over_jro_0_%i.png'% (time.time()) #plot pattern & objects | |||
|
1380 | self.plotname1 = 'over_jro_1_%i.png'% (time.time()) #plot antenna cuts | |||
|
1381 | self.plotname2 = 'over_jro_2_%i.png'% (time.time()) #plot sky noise | |||
|
1382 | ||||
|
1383 | # Defining antenna axes respect to geographic coordinates (See Ochs report). | |||
|
1384 | # alfa = 1.46*Misc_Routines.CoFactors.d2r | |||
|
1385 | # theta = 51.01*Misc_Routines.CoFactors.d2r | |||
|
1386 | ||||
|
1387 | alfa = 1.488312*Misc_Routines.CoFactors.d2r | |||
|
1388 | th = 6.166710 + 45.0 | |||
|
1389 | theta = th*Misc_Routines.CoFactors.d2r | |||
|
1390 | ||||
|
1391 | self.maxha_min = 4*self.maxphi*numpy.sqrt(2)*1.25 | |||
|
1392 | ||||
|
1393 | sina = numpy.sin(alfa) | |||
|
1394 | cosa = numpy.cos(alfa) | |||
|
1395 | MT1 = numpy.array([[1,0,0],[0,cosa,-sina],[0,sina,cosa]]) | |||
|
1396 | sinb = numpy.sin(theta) | |||
|
1397 | cosb = numpy.cos(theta) | |||
|
1398 | MT2 = numpy.array([[cosb,sinb,0],[-sinb,cosb,0],[0,0,1]]) | |||
|
1399 | self.MT3 = numpy.array(numpy.dot(MT2, MT1)).transpose() | |||
|
1400 | ||||
|
1401 | self.xg = numpy.dot(self.MT3.transpose(),numpy.array([1,0,0])) | |||
|
1402 | self.yg = numpy.dot(self.MT3.transpose(),numpy.array([0,1,0])) | |||
|
1403 | self.zg = numpy.dot(self.MT3.transpose(),numpy.array([0,0,1])) | |||
|
1404 | ||||
|
1405 | def plotPattern2(self, date, phases, gain_tx, gain_rx, ues, just_rx, angle): | |||
|
1406 | # Plotting Antenna patterns. | |||
|
1407 | self.maxphi = angle | |||
|
1408 | self.initParameters(date) | |||
|
1409 | self.doy = datetime.datetime(date.year,date.month,date.day).timetuple().tm_yday | |||
|
1410 | self.junkjd = TimeTools.Time(self.year,self.month,self.dom).change2julday() | |||
|
1411 | self.junklst = TimeTools.Julian(self.junkjd).change2lst(longitude=self.glon) | |||
|
1412 | self.ra_obs = self.junklst*Misc_Routines.CoFactors.h2d | |||
|
1413 | ||||
|
1414 | date = TimeTools.Time(date.year,date.month,date.day).change2strdate(mode=2) | |||
|
1415 | ||||
|
1416 | mesg = 'Over Jicamarca: ' + date[0] | |||
|
1417 | ||||
|
1418 | ObjAnt = JroPattern(pattern=0, | |||
|
1419 | filename=None, | |||
|
1420 | path=None, | |||
|
1421 | nptsx=self.nptsx, | |||
|
1422 | nptsy=self.nptsy, | |||
|
1423 | maxphi=angle, | |||
|
1424 | fftopt=self.fftopt, | |||
|
1425 | phases=phases, | |||
|
1426 | gain_tx=gain_tx, | |||
|
1427 | gain_rx=gain_rx, | |||
|
1428 | ues=ues, | |||
|
1429 | just_rx=just_rx | |||
|
1430 | ) | |||
|
1431 | ||||
|
1432 | self.pattern_plot = AntPatternPlot() | |||
|
1433 | ||||
|
1434 | self.pattern_plot.contPattern(iplot=0, | |||
|
1435 | gpath=self.path4plotname, | |||
|
1436 | filename=self.plotname0, | |||
|
1437 | mesg=mesg, | |||
|
1438 | amp=ObjAnt.norpattern, | |||
|
1439 | x=ObjAnt.dcosx, | |||
|
1440 | y=ObjAnt.dcosy, | |||
|
1441 | getCut=ObjAnt.getcut, | |||
|
1442 | title=self.ptitle, | |||
|
1443 | save=False) | |||
|
1444 | ||||
|
1445 | ||||
|
1446 | self.pattern_plot.plotRaDec(gpath=self.path4plotname, | |||
|
1447 | filename=self.plotname0, | |||
|
1448 | jd=self.junkjd, | |||
|
1449 | ra_obs=self.ra_obs, | |||
|
1450 | xg=self.xg, | |||
|
1451 | yg=self.yg, | |||
|
1452 | x=ObjAnt.dcosx, | |||
|
1453 | y=ObjAnt.dcosy, | |||
|
1454 | save=False) | |||
|
1455 | ||||
|
1456 | vect_ant = numpy.array([ObjAnt.meanpos[0],ObjAnt.meanpos[1],numpy.sqrt(1-numpy.sum(ObjAnt.meanpos**2.))]) | |||
|
1457 | vect_geo = numpy.dot(scipy.linalg.inv(self.MT3),vect_ant) | |||
|
1458 | vect_polar = Misc_Routines.Vector(numpy.array(vect_geo),direction=1).Polar2Rect() | |||
|
1459 | [ra,dec,ha] = Astro_Coords.AltAz(vect_polar[1],vect_polar[0],self.junkjd).change2equatorial() | |||
|
1460 | self.main_dec = dec | |||
|
1461 | ||||
|
1462 | def plotBfield(self): | |||
|
1463 | ||||
|
1464 | ObjB = BField(self.year,self.doy,1,self.heights) | |||
|
1465 | [dcos, alpha, nlon, nlat] = ObjB.getBField() | |||
|
1466 | ||||
|
1467 | self.pattern_plot.plotBField('', '',dcos,alpha,nlon,nlat, | |||
|
1468 | self.dcosxrange, | |||
|
1469 | self.dcosyrange, | |||
|
1470 | ObjB.heights, | |||
|
1471 | ObjB.alpha_i, | |||
|
1472 | save=False) | |||
|
1473 | ||||
|
1474 | def plotCelestial(self, objects): | |||
|
1475 | ||||
|
1476 | ntod = 24.*16. | |||
|
1477 | ||||
|
1478 | tod = numpy.arange(ntod)/ntod*24. | |||
|
1479 | ||||
|
1480 | [month,dom] = TimeTools.Doy2Date(self.year, self.doy).change2date() | |||
|
1481 | ||||
|
1482 | jd = TimeTools.Time(self.year, month, dom, tod+self.UT).change2julday() | |||
|
1483 | ||||
|
1484 | self.pattern_plot.plotCelestial( | |||
|
1485 | jd, | |||
|
1486 | self.main_dec, | |||
|
1487 | tod, | |||
|
1488 | self.maxha_min, | |||
|
1489 | objects, | |||
|
1490 | self.glat, | |||
|
1491 | self.glon, | |||
|
1492 | self.xg, | |||
|
1493 | self.yg, | |||
|
1494 | self.dcosxrange, | |||
|
1495 | self.dcosyrange | |||
|
1496 | ) | |||
|
1497 | ||||
|
1498 | def plotAntennaCuts(self): | |||
|
1499 | # print "Drawing antenna cuts" | |||
|
1500 | ||||
|
1501 | incha = 0.05 # min | |||
|
1502 | nha = numpy.int32(2*self.maxha_min/incha) + 1. | |||
|
1503 | newha = numpy.arange(nha)/nha*2.*self.maxha_min - self.maxha_min | |||
|
1504 | nha_star = numpy.int32(200./incha) | |||
|
1505 | star_ha = (numpy.arange(nha_star) - (nha_star/2))*nha_star | |||
|
1506 | ||||
|
1507 | #Init ObjCut for PatternCutPlot() | |||
|
1508 | view_objects = numpy.where(self.show_object>0) | |||
|
1509 | subplots = len(view_objects[0]) | |||
|
1510 | ObjCut = Graphics_OverJro.PatternCutPlot(subplots) | |||
|
1511 | ||||
|
1512 | for io in (numpy.arange(5)): | |||
|
1513 | if self.show_object[io]==2: | |||
|
1514 | if io==0: | |||
|
1515 | if self.dcosx_mag.size!=0: | |||
|
1516 | dcosx = self.dcosx_mag | |||
|
1517 | dcosy = self.dcosy_mag | |||
|
1518 | dcosz = 1 - numpy.sqrt(dcosx**2. + dcosy**2.) | |||
|
1519 | ||||
|
1520 | # Finding rotation of B respec to antenna coords. | |||
|
1521 | [mm,bb] = scipy.polyfit(dcosx,dcosy,1) | |||
|
1522 | alfa = 0.0 | |||
|
1523 | theta = -1.*numpy.arctan(mm) | |||
|
1524 | sina = numpy.sin(alfa); cosa = numpy.cos(alfa) | |||
|
1525 | MT1 = [[1,0,0],[0,cosa,-sina],[0,sina,cosa]] | |||
|
1526 | MT1 = numpy.array(MT1) | |||
|
1527 | sinb = numpy.sin(theta); cosb = numpy.cos(theta) | |||
|
1528 | MT2 = [[cosb,sinb,0],[-sinb,cosb,0],[0,0,1]] | |||
|
1529 | MT2 = numpy.array(MT2) | |||
|
1530 | MT3_mag = numpy.dot(MT2, MT1) | |||
|
1531 | MT3_mag = numpy.array(MT3_mag).transpose() | |||
|
1532 | # Getting dcos respec to B coords | |||
|
1533 | vector = numpy.array([dcosx,dcosy,dcosz]) | |||
|
1534 | nvector = numpy.dot(MT3_mag,vector) | |||
|
1535 | nvector = numpy.array(nvector).transpose() | |||
|
1536 | ||||
|
1537 | ## print 'Rotation (deg) %f'%(theta/Misc_Routines.CoFactors.d2r) | |||
|
1538 | ||||
|
1539 | yoffset = numpy.sum(nvector[:,1])/nvector[:,1].size | |||
|
1540 | # print 'Dcosyoffset %f'%(yoffset) | |||
|
1541 | ||||
|
1542 | ha = self.ha_mag*4. | |||
|
1543 | time = self.time_mag | |||
|
1544 | width_star = 0.1 # half width in minutes | |||
|
1545 | otitle = 'B Perp. cut' | |||
|
1546 | # else: | |||
|
1547 | # print "No B perp. over Observatory" | |||
|
1548 | # | |||
|
1549 | # | |||
|
1550 | elif io==1: | |||
|
1551 | if self.ObjC.dcosx_sun.size!=0: | |||
|
1552 | dcosx = self.ObjC.dcosx_sun | |||
|
1553 | dcosy = self.ObjC.dcosy_sun | |||
|
1554 | ha = self.ObjC.ha_sun*4.0 | |||
|
1555 | time = self.ObjC.time_sun | |||
|
1556 | width_star = 2. # half width in minutes | |||
|
1557 | otitle = 'Sun cut' | |||
|
1558 | # else: | |||
|
1559 | # print "Sun is not passing over Observatory" | |||
|
1560 | ||||
|
1561 | elif io==2: | |||
|
1562 | if self.ObjC.dcosx_moon.size!=0: | |||
|
1563 | dcosx = self.ObjC.dcosx_moon | |||
|
1564 | dcosy = self.ObjC.dcosy_moon | |||
|
1565 | ha = self.ObjC.ha_moon*4 | |||
|
1566 | time = self.ObjC.time_moon | |||
|
1567 | m_distance = 404114.6 # distance to the Earth in km | |||
|
1568 | m_diameter = 1734.4 # diameter in km. | |||
|
1569 | width_star = numpy.arctan(m_distance/m_diameter) | |||
|
1570 | width_star = width_star/2./Misc_Routines.CoFactors.d2r*4. | |||
|
1571 | otitle = 'Moon cut' | |||
|
1572 | # else: | |||
|
1573 | # print "Moon is not passing over Observatory" | |||
|
1574 | ||||
|
1575 | elif io==3: | |||
|
1576 | if self.ObjC.dcosx_hydra.size!=0: | |||
|
1577 | dcosx = self.ObjC.dcosx_hydra | |||
|
1578 | dcosy = self.ObjC.dcosy_hydra | |||
|
1579 | ha = self.ObjC.ha_hydra*4. | |||
|
1580 | time = self.ObjC.time_hydra | |||
|
1581 | width_star = 0.25 # half width in minutes | |||
|
1582 | otitle = 'Hydra cut' | |||
|
1583 | # else: | |||
|
1584 | # print "Hydra is not passing over Observatory" | |||
|
1585 | ||||
|
1586 | elif io==4: | |||
|
1587 | if self.ObjC.dcosx_galaxy.size!=0: | |||
|
1588 | dcosx = self.ObjC.dcosx_galaxy | |||
|
1589 | dcosy = self.ObjC.dcosy_galaxy | |||
|
1590 | ha = self.ObjC.ha_galaxy*4. | |||
|
1591 | time = self.ObjC.time_galaxy | |||
|
1592 | width_star = 25. # half width in minutes | |||
|
1593 | otitle = 'Galaxy cut' | |||
|
1594 | # else: | |||
|
1595 | # print "Galaxy center is not passing over Jicamarca" | |||
|
1596 | # | |||
|
1597 | # | |||
|
1598 | hour = numpy.int32(time) | |||
|
1599 | mins = numpy.int32((time - hour)*60.) | |||
|
1600 | secs = numpy.int32(((time - hour)*60. - mins)*60.) | |||
|
1601 | ||||
|
1602 | ObjT = TimeTools.Time(self.year,self.month,self.dom,hour,mins,secs) | |||
|
1603 | subtitle = ObjT.change2strdate() | |||
|
1604 | ||||
|
1605 | star_cut = numpy.exp(-(star_ha/width_star)**2./2.) | |||
|
1606 | ||||
|
1607 | pol = scipy.polyfit(ha,dcosx,3.) | |||
|
1608 | polx = numpy.poly1d(pol); newdcosx = polx(newha) | |||
|
1609 | pol = scipy.polyfit(ha,dcosy,3.) | |||
|
1610 | poly = numpy.poly1d(pol);newdcosy = poly(newha) | |||
|
1611 | ||||
|
1612 | patterns = [] | |||
|
1613 | for icut in numpy.arange(self.pattern.size): | |||
|
1614 | # Getting Antenna cut. | |||
|
1615 | Obj = JroPattern(dcosx=newdcosx, | |||
|
1616 | dcosy=newdcosy, | |||
|
1617 | getcut=1, | |||
|
1618 | pattern=self.pattern[icut], | |||
|
1619 | path=self.path, | |||
|
1620 | filename=self.filename[icut]) | |||
|
1621 | ||||
|
1622 | Obj.getPattern() | |||
|
1623 | ||||
|
1624 | patterns.append(Obj.pattern) | |||
|
1625 | ||||
|
1626 | ||||
|
1627 | ObjCut.drawCut(io, | |||
|
1628 | patterns, | |||
|
1629 | self.pattern.size, | |||
|
1630 | newha, | |||
|
1631 | otitle, | |||
|
1632 | subtitle, | |||
|
1633 | self.ptitle) | |||
|
1634 | ||||
|
1635 | ObjCut.saveFig(self.path4plotname,self.plotname1) | |||
|
1636 | ||||
|
1637 | ||||
12 |
|
1638 | |||
13 | def skyNoise(jd, ut=-5.0, longitude=-76.87, filename='/app/utils/galaxy.txt'): |
|
1639 | def skyNoise(jd, ut=-5.0, longitude=-76.87, filename='/app/utils/galaxy.txt'): | |
14 | """ |
|
1640 | """ | |
@@ -186,3 +1812,38 def skynoise_plot(year, month, day): | |||||
186 |
|
1812 | |||
187 | return buf |
|
1813 | return buf | |
188 |
|
1814 | |||
|
1815 | def overjro_plot(pattern_id, date, angle, height, bodys): | |||
|
1816 | ||||
|
1817 | pattern = select_pattern(pattern = pattern_id) | |||
|
1818 | ob = overJroShow(pattern['title'], heights=height) | |||
|
1819 | ob.plotPattern2( | |||
|
1820 | date, | |||
|
1821 | pattern['phase'], | |||
|
1822 | pattern['gaintx'], | |||
|
1823 | pattern['gainrx'], | |||
|
1824 | pattern['ues'], | |||
|
1825 | pattern['justrx'], | |||
|
1826 | angle | |||
|
1827 | ) | |||
|
1828 | ||||
|
1829 | if 'bfield' in bodys: | |||
|
1830 | ob.plotBfield() | |||
|
1831 | ||||
|
1832 | objects = [] | |||
|
1833 | for body in bodys: | |||
|
1834 | if body=='sun': | |||
|
1835 | objects.append(0) | |||
|
1836 | elif body=='moon': | |||
|
1837 | objects.append(1) | |||
|
1838 | elif body=='hydra': | |||
|
1839 | objects.append(2) | |||
|
1840 | elif body=='galaxy': | |||
|
1841 | objects.append(3) | |||
|
1842 | ||||
|
1843 | if objects: | |||
|
1844 | ob.plotCelestial(objects) | |||
|
1845 | ||||
|
1846 | buf = BytesIO() | |||
|
1847 | ob.pattern_plot.fig.savefig(buf, format="png") | |||
|
1848 | ||||
|
1849 | return buf No newline at end of file |
General Comments 0
You need to be logged in to leave comments.
Login now