Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1from __future__ import division, with_statement 

2from builtins import zip 

3from builtins import range 

4import numpy as np 

5import unittest 

6import lsst.utils.tests 

7import lsst.sims.utils as utils 

8 

9 

10def setup_module(module): 

11 lsst.utils.tests.init() 

12 

13 

14def controlEquationOfEquinoxes(mjd): 

15 """ 

16 Taken from http://aa.usno.navy.mil/faq/docs/GAST.php 

17 

18 @param [in] mjd is Terrestrial Time as a Modified Julian Date 

19 

20 @param [out] the equation of equinoxes in radians 

21 """ 

22 

23 JD = mjd + 2400000.5 

24 D = JD - 2451545.0 

25 omegaDegrees = 125.04 - 0.052954*D 

26 Ldegrees = 280.47 + 0.98565*D 

27 deltaPsiHours = -0.000319 * np.sin(np.radians(omegaDegrees)) \ 

28 - 0.000024 * np.sin(2.0*np.radians(Ldegrees)) 

29 epsilonDegrees = 23.4393 - 0.0000004*D 

30 return (deltaPsiHours / 24.0) * 2.0 * np.pi * np.cos(np.radians(epsilonDegrees)) 

31 

32 

33def controlCalcGmstGast(mjd): 

34 # From http://aa.usno.navy.mil/faq/docs/GAST.php Nov. 9 2013 

35 mjdConv = 2400000.5 

36 jd2000 = 2451545.0 

37 mjd_o = np.floor(mjd) 

38 jd = mjd + mjdConv 

39 jd_o = mjd_o + mjdConv 

40 h = 24. * (jd - jd_o) 

41 d = jd - jd2000 

42 d_o = jd_o - jd2000 

43 t = d / 36525. 

44 gmst = 6.697374558 + 0.06570982441908*d_o + 1.00273790935 * h + 0.000026 * t**2 

45 gast = gmst + 24.0 * utils.equationOfEquinoxes(mjd) / (2.0 * np.pi) 

46 gmst %= 24. 

47 gast %= 24. 

48 return gmst, gast 

49 

50 

51class AngularSeparationTestCase(unittest.TestCase): 

52 

53 def testAngSepExceptions(self): 

54 """ 

55 Test that an exception is raised when you pass 

56 mismatched inputs to angularSeparation. 

57 """ 

58 ra1 = 23.0 

59 dec1 = -12.0 

60 ra2 = 45.0 

61 dec2 = 33.1 

62 ra1_arr = np.array([11.0, 21.0, 33.1]) 

63 dec1_arr = np.array([-11.1, 34.1, 86.2]) 

64 ra2_arr = np.array([45.2, 112.0, 89.3]) 

65 dec2_arr = np.array([11.1, -45.0, -71.0]) 

66 

67 # test that everything runs 

68 utils._angularSeparation(np.radians(ra1), np.radians(dec1), 

69 np.radians(ra2), np.radians(dec2)) 

70 utils.angularSeparation(ra1, dec1, ra2, dec2) 

71 ans = utils._angularSeparation(np.radians(ra1_arr), np.radians(dec1_arr), 

72 np.radians(ra2_arr), np.radians(dec2_arr)) 

73 self.assertEqual(len(ans), 3) 

74 ans = utils.angularSeparation(ra1_arr, dec1_arr, ra2_arr, dec2_arr) 

75 self.assertEqual(len(ans), 3) 

76 

77 ans = utils.angularSeparation(ra1_arr, dec1_arr, ra2, dec2) 

78 self.assertEqual(len(ans), 3) 

79 ans = utils._angularSeparation(ra1_arr, ra2_arr, ra2, dec2) 

80 self.assertEqual(len(ans), 3) 

81 

82 ans = utils.angularSeparation(ra1, dec1, ra2_arr, dec2_arr) 

83 self.assertEqual(len(ans), 3) 

84 ans = utils._angularSeparation(dec1, ra1, ra2_arr, dec2_arr) 

85 self.assertEqual(len(ans), 3) 

86 

87 # test with lists 

88 # Note: these exceptions will not get raised if you call 

89 # angularSeparation() with lists instead of numpy arrays 

90 # because the conversion from degrees to radians with 

91 # np.radians will automatically convert the lists into arrays 

92 with self.assertRaises(RuntimeError) as context: 

93 utils._angularSeparation(list(ra1_arr), dec1_arr, ra2_arr, dec2_arr) 

94 self.assertIn("number", context.exception.args[0]) 

95 self.assertIn("numpy array", context.exception.args[0]) 

96 self.assertIn("long1", context.exception.args[0]) 

97 self.assertIn("angularSeparation", context.exception.args[0]) 

98 

99 with self.assertRaises(RuntimeError) as context: 

100 utils._angularSeparation(ra1_arr, list(dec1_arr), ra2_arr, dec2_arr) 

101 self.assertIn("number", context.exception.args[0]) 

102 self.assertIn("numpy array", context.exception.args[0]) 

103 self.assertIn("lat1", context.exception.args[0]) 

104 self.assertIn("angularSeparation", context.exception.args[0]) 

105 

106 with self.assertRaises(RuntimeError) as context: 

107 utils._angularSeparation(ra1_arr, dec1_arr, list(ra2_arr), dec2_arr) 

108 self.assertIn("number", context.exception.args[0]) 

109 self.assertIn("numpy array", context.exception.args[0]) 

110 self.assertIn("long2", context.exception.args[0]) 

111 self.assertIn("angularSeparation", context.exception.args[0]) 

112 

113 with self.assertRaises(RuntimeError) as context: 

114 utils._angularSeparation(ra1_arr, dec1_arr, ra2_arr, list(dec2_arr)) 

115 self.assertIn("number", context.exception.args[0]) 

116 self.assertIn("numpy array", context.exception.args[0]) 

117 self.assertIn("lat2", context.exception.args[0]) 

118 self.assertIn("angularSeparation", context.exception.args[0]) 

119 

120 # test with numbers and arrays 

121 with self.assertRaises(RuntimeError) as context: 

122 utils._angularSeparation(ra1_arr, dec1, ra2, dec2) 

123 self.assertIn("the same type", context.exception.args[0]) 

124 self.assertIn("angularSeparation", context.exception.args[0]) 

125 

126 with self.assertRaises(RuntimeError) as context: 

127 utils._angularSeparation(ra1, dec1_arr, ra2, dec2) 

128 self.assertIn("the same type", context.exception.args[0]) 

129 self.assertIn("angularSeparation", context.exception.args[0]) 

130 

131 with self.assertRaises(RuntimeError) as context: 

132 utils._angularSeparation(ra1, dec1, ra2_arr, dec2) 

133 self.assertIn("the same type", context.exception.args[0]) 

134 self.assertIn("angularSeparation", context.exception.args[0]) 

135 

136 with self.assertRaises(RuntimeError) as context: 

137 utils._angularSeparation(ra1, dec1, ra2, dec2_arr) 

138 self.assertIn("the same type", context.exception.args[0]) 

139 self.assertIn("angularSeparation", context.exception.args[0]) 

140 

141 with self.assertRaises(RuntimeError) as context: 

142 utils.angularSeparation(ra1_arr, dec1, ra2, dec2) 

143 self.assertIn("the same type", context.exception.args[0]) 

144 self.assertIn("angularSeparation", context.exception.args[0]) 

145 

146 with self.assertRaises(RuntimeError) as context: 

147 utils.angularSeparation(ra1, dec1_arr, ra2, dec2) 

148 self.assertIn("the same type", context.exception.args[0]) 

149 self.assertIn("angularSeparation", context.exception.args[0]) 

150 

151 with self.assertRaises(RuntimeError) as context: 

152 utils.angularSeparation(ra1, dec1, ra2_arr, dec2) 

153 self.assertIn("the same type", context.exception.args[0]) 

154 self.assertIn("angularSeparation", context.exception.args[0]) 

155 

156 with self.assertRaises(RuntimeError) as context: 

157 utils.angularSeparation(ra1, dec1, ra2, dec2_arr) 

158 self.assertIn("the same type", context.exception.args[0]) 

159 self.assertIn("angularSeparation", context.exception.args[0]) 

160 

161 # test with mismatched arrays 

162 with self.assertRaises(RuntimeError) as context: 

163 utils._angularSeparation(ra1_arr[:2], dec1_arr, ra2_arr, dec2_arr) 

164 self.assertIn("same length", context.exception.args[0]) 

165 self.assertIn("angularSeparation", context.exception.args[0]) 

166 

167 with self.assertRaises(RuntimeError) as context: 

168 utils._angularSeparation(ra1_arr, dec1_arr[:2], ra2_arr, dec2_arr) 

169 self.assertIn("same length", context.exception.args[0]) 

170 self.assertIn("angularSeparation", context.exception.args[0]) 

171 

172 with self.assertRaises(RuntimeError) as context: 

173 utils._angularSeparation(ra1_arr, dec1_arr, ra2_arr[:2], dec2_arr) 

174 self.assertIn("same length", context.exception.args[0]) 

175 self.assertIn("angularSeparation", context.exception.args[0]) 

176 

177 with self.assertRaises(RuntimeError) as context: 

178 utils._angularSeparation(ra1_arr, dec1_arr, ra2_arr, dec2_arr[:2]) 

179 self.assertIn("same length", context.exception.args[0]) 

180 self.assertIn("angularSeparation", context.exception.args[0]) 

181 

182 with self.assertRaises(RuntimeError) as context: 

183 utils.angularSeparation(ra1_arr[:2], dec1_arr, ra2_arr, dec2_arr) 

184 self.assertIn("same length", context.exception.args[0]) 

185 self.assertIn("angularSeparation", context.exception.args[0]) 

186 

187 with self.assertRaises(RuntimeError) as context: 

188 utils.angularSeparation(ra1_arr, dec1_arr[:2], ra2_arr, dec2_arr) 

189 self.assertIn("same length", context.exception.args[0]) 

190 self.assertIn("angularSeparation", context.exception.args[0]) 

191 

192 with self.assertRaises(RuntimeError) as context: 

193 utils.angularSeparation(ra1_arr, dec1_arr, ra2_arr[:2], dec2_arr) 

194 self.assertIn("same length", context.exception.args[0]) 

195 self.assertIn("angularSeparation", context.exception.args[0]) 

196 

197 with self.assertRaises(RuntimeError) as context: 

198 utils.angularSeparation(ra1_arr, dec1_arr, ra2_arr, dec2_arr[:2]) 

199 self.assertIn("same length", context.exception.args[0]) 

200 self.assertIn("angularSeparation", context.exception.args[0]) 

201 

202 with self.assertRaises(RuntimeError) as context: 

203 utils._angularSeparation(ra1_arr[:2], dec1_arr[:2], ra2_arr, dec2_arr) 

204 self.assertIn("same length", context.exception.args[0]) 

205 self.assertIn("angularSeparation", context.exception.args[0]) 

206 

207 with self.assertRaises(RuntimeError) as context: 

208 utils._angularSeparation(ra1_arr, dec1_arr, ra2_arr[:2], dec2_arr[:2]) 

209 self.assertIn("same length", context.exception.args[0]) 

210 self.assertIn("angularSeparation", context.exception.args[0]) 

211 

212 with self.assertRaises(RuntimeError) as context: 

213 utils.angularSeparation(ra1_arr[:2], dec1_arr[:2], ra2_arr, dec2_arr) 

214 self.assertIn("same length", context.exception.args[0]) 

215 self.assertIn("angularSeparation", context.exception.args[0]) 

216 

217 with self.assertRaises(RuntimeError) as context: 

218 utils.angularSeparation(ra1_arr, dec1_arr, ra2_arr[:2], dec2_arr[:2]) 

219 self.assertIn("same length", context.exception.args[0]) 

220 self.assertIn("angularSeparation", context.exception.args[0]) 

221 

222 # test that a sensible error is raised if you pass a string 

223 # into angularSeparation 

224 # Note: a different error will be raised if you pass these 

225 # bad inputs into angularSeparation(). The exception will come 

226 # from trying to convert a str with np.radians 

227 with self.assertRaises(RuntimeError) as context: 

228 utils._angularSeparation('a', dec1, ra2, dec2) 

229 self.assertIn("angularSeparation", context.exception.args[0]) 

230 self.assertIn("number", context.exception.args[0]) 

231 

232 with self.assertRaises(RuntimeError) as context: 

233 utils._angularSeparation(ra1, 'a', ra2, dec2) 

234 self.assertIn("angularSeparation", context.exception.args[0]) 

235 self.assertIn("number", context.exception.args[0]) 

236 

237 with self.assertRaises(RuntimeError) as context: 

238 utils._angularSeparation(ra1, dec1, 'a', dec2) 

239 self.assertIn("angularSeparation", context.exception.args[0]) 

240 self.assertIn("number", context.exception.args[0]) 

241 

242 with self.assertRaises(RuntimeError) as context: 

243 utils._angularSeparation(ra1, dec1, ra2, 'a') 

244 self.assertIn("angularSeparation", context.exception.args[0]) 

245 self.assertIn("number", context.exception.args[0]) 

246 

247 def testAngSepResultsArr(self): 

248 """ 

249 Test that angularSeparation gives the correct answer by comparing 

250 results with the dot products of Cartesian vectors. Pass in arrays 

251 of arguments. 

252 """ 

253 rng = np.random.RandomState(99421) 

254 n_obj = 100 

255 ra1 = rng.random_sample(n_obj)*2.0*np.pi 

256 dec1 = rng.random_sample(n_obj)*np.pi - 0.5*np.pi 

257 ra2 = rng.random_sample(n_obj)*2.0*np.pi 

258 dec2 = rng.random_sample(n_obj)*np.pi - 0.5*np.pi 

259 

260 x1 = np.cos(dec1)*np.cos(ra1) 

261 y1 = np.cos(dec1)*np.sin(ra1) 

262 z1 = np.sin(dec1) 

263 

264 x2 = np.cos(dec2)*np.cos(ra2) 

265 y2 = np.cos(dec2)*np.sin(ra2) 

266 z2 = np.sin(dec2) 

267 

268 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

269 test = np.cos(test) 

270 control = x1*x2 + y1*y2 + z1*z2 

271 np.testing.assert_array_almost_equal(test, control, decimal=15) 

272 

273 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

274 np.degrees(ra2), np.degrees(dec2)) 

275 test = np.cos(np.radians(test)) 

276 np.testing.assert_array_almost_equal(test, control, decimal=15) 

277 

278 # specifically test at the north pole 

279 dec1 = np.ones(n_obj)*np.pi 

280 x1 = np.cos(dec1)*np.cos(ra1) 

281 y1 = np.cos(dec1)*np.sin(ra1) 

282 z1 = np.sin(dec1) 

283 control = x1*x2 + y1*y2 + z1*z2 

284 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

285 test = np.cos(test) 

286 np.testing.assert_array_almost_equal(test, control, decimal=15) 

287 

288 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

289 np.degrees(ra2), np.degrees(dec2)) 

290 test = np.cos(np.radians(test)) 

291 np.testing.assert_array_almost_equal(test, control, decimal=15) 

292 

293 # specifically test at the south pole 

294 dec1 = -1.0*np.ones(n_obj)*np.pi 

295 x1 = np.cos(dec1)*np.cos(ra1) 

296 y1 = np.cos(dec1)*np.sin(ra1) 

297 z1 = np.sin(dec1) 

298 control = x1*x2 + y1*y2 + z1*z2 

299 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

300 test = np.cos(test) 

301 np.testing.assert_array_almost_equal(test, control, decimal=15) 

302 

303 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

304 np.degrees(ra2), np.degrees(dec2)) 

305 test = np.cos(np.radians(test)) 

306 np.testing.assert_array_almost_equal(test, control, decimal=15) 

307 

308 def testAngSepResultsExtreme(self): 

309 """ 

310 Test that angularSeparation gives the correct answer by comparing 

311 results with the dot products of Cartesian vectors. Test on extremal 

312 values (i.e. longitudes that go beyond 360.0 and latitudes that go 

313 beyond 90.0) 

314 """ 

315 rng = np.random.RandomState(99421) 

316 n_obj = 100 

317 for sgn in (-1.0, 1.0): 

318 ra1 = sgn*(rng.random_sample(n_obj)*2.0*np.pi + 2.0*np.pi) 

319 dec1 = sgn*(rng.random_sample(n_obj)*4.0*np.pi + 2.0*np.pi) 

320 ra2 = sgn*(rng.random_sample(n_obj)*2.0*np.pi + 2.0*np.pi) 

321 dec2 = sgn*(rng.random_sample(n_obj)*2.0*np.pi + 2.0*np.pi) 

322 

323 x1 = np.cos(dec1)*np.cos(ra1) 

324 y1 = np.cos(dec1)*np.sin(ra1) 

325 z1 = np.sin(dec1) 

326 

327 x2 = np.cos(dec2)*np.cos(ra2) 

328 y2 = np.cos(dec2)*np.sin(ra2) 

329 z2 = np.sin(dec2) 

330 

331 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

332 test = np.cos(test) 

333 control = x1*x2 + y1*y2 + z1*z2 

334 np.testing.assert_array_almost_equal(test, control, decimal=14) 

335 

336 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

337 np.degrees(ra2), np.degrees(dec2)) 

338 test = np.cos(np.radians(test)) 

339 np.testing.assert_array_almost_equal(test, control, decimal=14) 

340 

341 # specifically test at the north pole 

342 dec1 = np.ones(n_obj)*np.pi 

343 x1 = np.cos(dec1)*np.cos(ra1) 

344 y1 = np.cos(dec1)*np.sin(ra1) 

345 z1 = np.sin(dec1) 

346 control = x1*x2 + y1*y2 + z1*z2 

347 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

348 test = np.cos(test) 

349 np.testing.assert_array_almost_equal(test, control, decimal=14) 

350 

351 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

352 np.degrees(ra2), np.degrees(dec2)) 

353 test = np.cos(np.radians(test)) 

354 dd = np.abs(test-control) 

355 np.testing.assert_array_almost_equal(test, control, decimal=14) 

356 

357 # specifically test at the south pole 

358 dec1 = -1.0*np.ones(n_obj)*np.pi 

359 x1 = np.cos(dec1)*np.cos(ra1) 

360 y1 = np.cos(dec1)*np.sin(ra1) 

361 z1 = np.sin(dec1) 

362 control = x1*x2 + y1*y2 + z1*z2 

363 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

364 test = np.cos(test) 

365 np.testing.assert_array_almost_equal(test, control, decimal=14) 

366 

367 test = utils.angularSeparation(np.degrees(ra1), np.degrees(dec1), 

368 np.degrees(ra2), np.degrees(dec2)) 

369 test = np.cos(np.radians(test)) 

370 np.testing.assert_array_almost_equal(test, control, decimal=14) 

371 

372 def testAngSepResultsFloat(self): 

373 """ 

374 Test that angularSeparation gives the correct answer by comparing 

375 results with the dot products of Cartesian vectors. Pass in floats 

376 as arguments. 

377 """ 

378 rng = np.random.RandomState(831) 

379 ra1 = rng.random_sample()*2.0*np.pi 

380 dec1 = rng.random_sample()*np.pi-0.5*np.pi 

381 ra2 = rng.random_sample()*2.0*np.pi 

382 dec2 = rng.random_sample()*np.pi-0.5*np.pi 

383 x1 = np.cos(dec1)*np.cos(ra1) 

384 y1 = np.cos(dec1)*np.sin(ra1) 

385 z1 = np.sin(dec1) 

386 

387 x2 = np.cos(dec2)*np.cos(ra2) 

388 y2 = np.cos(dec2)*np.sin(ra2) 

389 z2 = np.sin(dec2) 

390 

391 control = x1*x2 + y1*y2 + z1*z2 

392 

393 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

394 self.assertIsInstance(test, float) 

395 test = np.cos(test) 

396 self.assertAlmostEqual(control, test, 15) 

397 

398 test = utils._angularSeparation(np.array([ra1]), np.array([dec1]), ra2, dec2) 

399 self.assertIsInstance(test, float) 

400 test = np.cos(test) 

401 self.assertAlmostEqual(control, test, 15) 

402 

403 test = utils._angularSeparation(ra1, dec1, np.array([ra2]), np.array([dec2])) 

404 self.assertIsInstance(test, float) 

405 test = np.cos(test) 

406 self.assertAlmostEqual(control, test, 15) 

407 

408 # try north pole 

409 ra1 = 0.5*np.pi 

410 x1 = np.cos(dec1)*np.cos(ra1) 

411 y1 = np.cos(dec1)*np.sin(ra1) 

412 z1 = np.sin(dec1) 

413 control = x1*x2 + y1*y2 + z1*z2 

414 

415 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

416 self.assertIsInstance(test, float) 

417 test = np.cos(test) 

418 self.assertAlmostEqual(control, test, 15) 

419 

420 test = utils._angularSeparation(np.array([ra1]), np.array([dec1]), ra2, dec2) 

421 self.assertIsInstance(test, float) 

422 test = np.cos(test) 

423 self.assertAlmostEqual(control, test, 15) 

424 

425 test = utils._angularSeparation(ra1, dec1, np.array([ra2]), np.array([dec2])) 

426 self.assertIsInstance(test, float) 

427 test = np.cos(test) 

428 self.assertAlmostEqual(control, test, 15) 

429 

430 # do all of that in degrees 

431 ra1 = rng.random_sample()*360.0 

432 dec1 = rng.random_sample()*180.0-90.0 

433 ra2 = rng.random_sample()*360.0 

434 dec2 = rng.random_sample()*180.0-90.0 

435 x1 = np.cos(np.radians(dec1))*np.cos(np.radians(ra1)) 

436 y1 = np.cos(np.radians(dec1))*np.sin(np.radians(ra1)) 

437 z1 = np.sin(np.radians(dec1)) 

438 

439 x2 = np.cos(np.radians(dec2))*np.cos(np.radians(ra2)) 

440 y2 = np.cos(np.radians(dec2))*np.sin(np.radians(ra2)) 

441 z2 = np.sin(np.radians(dec2)) 

442 

443 control = x1*x2 + y1*y2 + z1*z2 

444 

445 test = utils.angularSeparation(ra1, dec1, ra2, dec2) 

446 self.assertIsInstance(test, float) 

447 test = np.cos(np.radians(test)) 

448 self.assertAlmostEqual(control, test, 15) 

449 

450 test = utils.angularSeparation(np.array([ra1]), np.array([dec1]), ra2, dec2) 

451 self.assertIsInstance(test, float) 

452 test = np.cos(np.radians(test)) 

453 self.assertAlmostEqual(control, test, 15) 

454 

455 test = utils.angularSeparation(ra1, dec1, np.array([ra2]), np.array([dec2])) 

456 self.assertIsInstance(test, float) 

457 test = np.cos(np.radians(test)) 

458 self.assertAlmostEqual(control, test, 15) 

459 

460 # try north pole 

461 ra1 = 90.0 

462 x1 = np.cos(np.radians(dec1))*np.cos(np.radians(ra1)) 

463 y1 = np.cos(np.radians(dec1))*np.sin(np.radians(ra1)) 

464 z1 = np.sin(np.radians(dec1)) 

465 control = x1*x2 + y1*y2 + z1*z2 

466 

467 test = utils.angularSeparation(ra1, dec1, ra2, dec2) 

468 self.assertIsInstance(test, float) 

469 test = np.cos(np.radians(test)) 

470 self.assertAlmostEqual(control, test, 15) 

471 

472 test = utils.angularSeparation(np.array([ra1]), np.array([dec1]), ra2, dec2) 

473 self.assertIsInstance(test, float) 

474 test = np.cos(np.radians(test)) 

475 self.assertAlmostEqual(control, test, 15) 

476 

477 test = utils.angularSeparation(ra1, dec1, np.array([ra2]), np.array([dec2])) 

478 self.assertIsInstance(test, float) 

479 test = np.cos(np.radians(test)) 

480 self.assertAlmostEqual(control, test, 15) 

481 

482 def testAngSepResultsMixed(self): 

483 """ 

484 Test that angularSeparation gives the correct answer by comparing 

485 results with the dot products of Cartesian vectors. Pass in mixtures 

486 of floats and arrays as arguments. 

487 """ 

488 rng = np.random.RandomState(8131) 

489 n_obj = 100 

490 ra1 = rng.random_sample(n_obj)*2.0*np.pi 

491 dec1 = rng.random_sample(n_obj)*np.pi - 0.5*np.pi 

492 ra2 = rng.random_sample()*2.0*np.pi 

493 dec2 = rng.random_sample()*np.pi - 0.5*np.pi 

494 self.assertIsInstance(ra1, np.ndarray) 

495 self.assertIsInstance(dec1, np.ndarray) 

496 self.assertIsInstance(ra2, float) 

497 self.assertIsInstance(dec2, float) 

498 

499 x1 = np.cos(dec1)*np.cos(ra1) 

500 y1 = np.cos(dec1)*np.sin(ra1) 

501 z1 = np.sin(dec1) 

502 

503 x2 = np.cos(dec2)*np.cos(ra2) 

504 y2 = np.cos(dec2)*np.sin(ra2) 

505 z2 = np.sin(dec2) 

506 

507 control = x1*x2 + y1*y2 + z1*z2 

508 test = utils._angularSeparation(ra1, dec1, ra2, dec2) 

509 test = np.cos(test) 

510 np.testing.assert_array_almost_equal(test, control, decimal=15) 

511 test = utils._angularSeparation(ra2, dec2, ra1, dec1) 

512 test = np.cos(test) 

513 np.testing.assert_array_almost_equal(test, control, decimal=15) 

514 

515 # now do it in degrees 

516 ra1 = rng.random_sample(n_obj)*360.0 

517 dec1 = rng.random_sample(n_obj)*180.0 - 90.0 

518 ra2 = rng.random_sample()*360.0 

519 dec2 = rng.random_sample()*180.0 - 90.0 

520 self.assertIsInstance(ra1, np.ndarray) 

521 self.assertIsInstance(dec1, np.ndarray) 

522 self.assertIsInstance(ra2, float) 

523 self.assertIsInstance(dec2, float) 

524 

525 x1 = np.cos(np.radians(dec1))*np.cos(np.radians(ra1)) 

526 y1 = np.cos(np.radians(dec1))*np.sin(np.radians(ra1)) 

527 z1 = np.sin(np.radians(dec1)) 

528 

529 x2 = np.cos(np.radians(dec2))*np.cos(np.radians(ra2)) 

530 y2 = np.cos(np.radians(dec2))*np.sin(np.radians(ra2)) 

531 z2 = np.sin(np.radians(dec2)) 

532 

533 control = x1*x2 + y1*y2 + z1*z2 

534 test = utils.angularSeparation(ra1, dec1, ra2, dec2) 

535 test = np.cos(np.radians(test)) 

536 np.testing.assert_array_almost_equal(test, control, decimal=15) 

537 test = utils.angularSeparation(ra2, dec2, ra1, dec1) 

538 test = np.cos(np.radians(test)) 

539 np.testing.assert_array_almost_equal(test, control, decimal=15) 

540 

541 def testHaversine(self): 

542 """ 

543 Test that haversine() returns the same thing as _angularSeparation 

544 """ 

545 

546 ra1 = 0.2 

547 dec1 = 1.3 

548 ra2 = 2.1 

549 dec2 = -0.5 

550 ra3 = np.array([1.9, 2.1, 0.3]) 

551 dec3 = np.array([-1.1, 0.34, 0.01]) 

552 control = utils._angularSeparation(ra1, dec1, ra2, dec2) 

553 test = utils.haversine(ra1, dec1, ra2, dec2) 

554 self.assertIsInstance(test, float) 

555 self.assertEqual(test, control) 

556 

557 control = utils._angularSeparation(ra1, dec1, ra3, dec3) 

558 test = utils.haversine(ra1, dec1, ra3, dec3) 

559 np.testing.assert_array_equal(test, control) 

560 

561 control = utils._angularSeparation(np.array([ra1]), np.array([dec1]), ra3, dec3) 

562 test = utils.haversine(np.array([ra1]), np.array([dec1]), ra3, dec3) 

563 np.testing.assert_array_equal(test, control) 

564 

565 control = utils._angularSeparation(ra3, dec3, np.array([ra1]), np.array([dec1])) 

566 test = utils.haversine(ra3, dec3, np.array([ra1]), np.array([dec1])) 

567 np.testing.assert_array_equal(test, control) 

568 

569 control = utils._angularSeparation(ra2, dec2, np.array([ra1]), np.array([dec1])) 

570 test = utils.haversine(ra2, dec2, np.array([ra1]), np.array([dec1])) 

571 self.assertIsInstance(test, float) 

572 self.assertEqual(test, control) 

573 

574 control = utils._angularSeparation(np.array([ra1]), np.array([dec1]), ra2, dec2) 

575 test = utils.haversine(np.array([ra1]), np.array([dec1]), ra2, dec2) 

576 self.assertIsInstance(test, float) 

577 self.assertEqual(test, control) 

578 

579 

580class testCoordinateTransformations(unittest.TestCase): 

581 

582 def setUp(self): 

583 self.rng = np.random.RandomState(32) 

584 ntests = 100 

585 self.mjd = 57087.0 - 1000.0 * (self.rng.random_sample(ntests) - 0.5) 

586 self.tolerance = 1.0e-5 

587 

588 def testExceptions(self): 

589 """ 

590 Test to make sure that methods complain when incorrect data types are passed. 

591 """ 

592 mjdFloat = 52000.0 

593 mjd2 = np.array([52000.0, 53000.0]) 

594 mjd3 = np.array([53000.0, 53000.0, 54000.0]) 

595 

596 longFloat = 1.2 

597 longArr = np.array([1.2, 1.4]) 

598 

599 self.assertRaises(RuntimeError, utils.calcLmstLast, mjdFloat, longArr) 

600 self.assertRaises(RuntimeError, utils.calcLmstLast, mjd3, longArr) 

601 self.assertRaises(RuntimeError, utils.calcLmstLast, list(mjd2), longArr) 

602 self.assertRaises(RuntimeError, utils.calcLmstLast, mjd2, list(longArr)) 

603 self.assertRaises(RuntimeError, utils.calcLmstLast, mjdFloat, longArr) 

604 utils.calcLmstLast(mjd2, longFloat) 

605 utils.calcLmstLast(mjdFloat, longFloat) 

606 utils.calcLmstLast(int(mjdFloat), longFloat) 

607 utils.calcLmstLast(mjdFloat, int(longFloat)) 

608 utils.calcLmstLast(int(mjdFloat), int(longFloat)) 

609 utils.calcLmstLast(mjd2, longArr) 

610 

611 def testEquationOfEquinoxes(self): 

612 """ 

613 Test equation of equninoxes calculation 

614 """ 

615 

616 # test vectorized version 

617 control = controlEquationOfEquinoxes(self.mjd) 

618 test = utils.equationOfEquinoxes(self.mjd) 

619 self.assertLess(np.abs(test-control).max(), self.tolerance) 

620 

621 # test non-vectorized version 

622 for mm in self.mjd: 

623 control = controlEquationOfEquinoxes(mm) 

624 test = utils.equationOfEquinoxes(mm) 

625 self.assertLess(np.abs(test-control), self.tolerance) 

626 

627 def testGmstGast(self): 

628 """ 

629 Test calculation of Greenwich mean and apparent sidereal times 

630 """ 

631 

632 controlGmst, controlGast = controlCalcGmstGast(self.mjd) 

633 testGmst, testGast = utils.calcGmstGast(self.mjd) 

634 self.assertLess(np.abs(testGmst - controlGmst).max(), self.tolerance) 

635 self.assertLess(np.abs(testGast - controlGast).max(), self.tolerance) 

636 

637 # test non-vectorized version 

638 for mm in self.mjd: 

639 controlGmst, controlGast = controlCalcGmstGast(mm) 

640 testGmst, testGast = utils.calcGmstGast(mm) 

641 self.assertLess(np.abs(testGmst - controlGmst), self.tolerance) 

642 self.assertLess(np.abs(testGast - controlGast), self.tolerance) 

643 

644 def testLmstLast(self): 

645 """ 

646 Test calculation of local mean and apparent sidereal time 

647 """ 

648 

649 gmst, gast = utils.calcGmstGast(self.mjd) 

650 ll = [1.2, 2.2] 

651 

652 # test passing a float for longitude and a numpy array for mjd 

653 for longitude in ll: 

654 hours = np.degrees(longitude) / 15.0 

655 if hours > 24.0: 

656 hours -= 24.0 

657 controlLmst = gmst + hours 

658 controlLast = gast + hours 

659 controlLmst %= 24.0 

660 controlLast %= 24.0 

661 testLmst, testLast = utils.calcLmstLast(self.mjd, longitude) 

662 self.assertLess(np.abs(testLmst - controlLmst).max(), self.tolerance) 

663 self.assertLess(np.abs(testLast - controlLast).max(), self.tolerance) 

664 self.assertIsInstance(testLmst, np.ndarray) 

665 self.assertIsInstance(testLast, np.ndarray) 

666 

667 # test passing two floats 

668 for longitude in ll: 

669 for mm in self.mjd: 

670 gmst, gast = utils.calcGmstGast(mm) 

671 hours = np.degrees(longitude) / 15.0 

672 if hours > 24.0: 

673 hours -= 24.0 

674 controlLmst = gmst + hours 

675 controlLast = gast + hours 

676 controlLmst %= 24.0 

677 controlLast %= 24.0 

678 testLmst, testLast = utils.calcLmstLast(mm, longitude) 

679 self.assertLess(np.abs(testLmst - controlLmst), self.tolerance) 

680 self.assertLess(np.abs(testLast - controlLast), self.tolerance) 

681 self.assertIsInstance(testLmst, np.float) 

682 self.assertIsInstance(testLast, np.float) 

683 

684 # test passing two numpy arrays 

685 ll = self.rng.random_sample(len(self.mjd)) * 2.0 * np.pi 

686 testLmst, testLast = utils.calcLmstLast(self.mjd, ll) 

687 self.assertIsInstance(testLmst, np.ndarray) 

688 self.assertIsInstance(testLast, np.ndarray) 

689 for ix, (longitude, mm) in enumerate(zip(ll, self.mjd)): 

690 controlLmst, controlLast = utils.calcLmstLast(mm, longitude) 

691 self.assertAlmostEqual(controlLmst, testLmst[ix], 10) 

692 self.assertAlmostEqual(controlLast, testLast[ix], 10) 

693 

694 def test_galacticFromEquatorial(self): 

695 

696 ra = np.zeros((3), dtype=float) 

697 dec = np.zeros((3), dtype=float) 

698 

699 ra[0] = 2.549091039839124218e+00 

700 dec[0] = 5.198752733024248895e-01 

701 ra[1] = 8.693375673649429425e-01 

702 dec[1] = 1.038086165642298164e+00 

703 ra[2] = 7.740864769302191473e-01 

704 dec[2] = 2.758053025017753179e-01 

705 

706 glon, glat = utils._galacticFromEquatorial(ra, dec) 

707 

708 self.assertIsInstance(glon, np.ndarray) 

709 self.assertIsInstance(glat, np.ndarray) 

710 

711 self.assertAlmostEqual(glon[0], 3.452036693523627964e+00, 6) 

712 self.assertAlmostEqual(glat[0], 8.559512505657201897e-01, 6) 

713 self.assertAlmostEqual(glon[1], 2.455968474619387720e+00, 6) 

714 self.assertAlmostEqual(glat[1], 3.158563770667878468e-02, 6) 

715 self.assertAlmostEqual(glon[2], 2.829585540991265358e+00, 6) 

716 self.assertAlmostEqual(glat[2], -6.510790587552289788e-01, 6) 

717 

718 # test passing in floats as args 

719 for ix, (rr, dd) in enumerate(zip(ra, dec)): 

720 gl, gb = utils._galacticFromEquatorial(rr, dd) 

721 self.assertIsInstance(rr, np.float) 

722 self.assertIsInstance(dd, np.float) 

723 self.assertIsInstance(gl, np.float) 

724 self.assertIsInstance(gb, np.float) 

725 self.assertAlmostEqual(gl, glon[ix], 10) 

726 self.assertAlmostEqual(gb, glat[ix], 10) 

727 

728 def test_equatorialFromGalactic(self): 

729 

730 lon = np.zeros((3), dtype=float) 

731 lat = np.zeros((3), dtype=float) 

732 

733 lon[0] = 3.452036693523627964e+00 

734 lat[0] = 8.559512505657201897e-01 

735 lon[1] = 2.455968474619387720e+00 

736 lat[1] = 3.158563770667878468e-02 

737 lon[2] = 2.829585540991265358e+00 

738 lat[2] = -6.510790587552289788e-01 

739 

740 ra, dec = utils._equatorialFromGalactic(lon, lat) 

741 

742 self.assertIsInstance(ra, np.ndarray) 

743 self.assertIsInstance(dec, np.ndarray) 

744 

745 self.assertAlmostEqual(ra[0], 2.549091039839124218e+00, 6) 

746 self.assertAlmostEqual(dec[0], 5.198752733024248895e-01, 6) 

747 self.assertAlmostEqual(ra[1], 8.693375673649429425e-01, 6) 

748 self.assertAlmostEqual(dec[1], 1.038086165642298164e+00, 6) 

749 self.assertAlmostEqual(ra[2], 7.740864769302191473e-01, 6) 

750 self.assertAlmostEqual(dec[2], 2.758053025017753179e-01, 6) 

751 

752 # test passing in floats as args 

753 for ix, (ll, bb) in enumerate(zip(lon, lat)): 

754 rr, dd = utils._equatorialFromGalactic(ll, bb) 

755 self.assertIsInstance(ll, np.float) 

756 self.assertIsInstance(bb, np.float) 

757 self.assertIsInstance(rr, np.float) 

758 self.assertIsInstance(dd, np.float) 

759 self.assertAlmostEqual(rr, ra[ix], 10) 

760 self.assertAlmostEqual(dd, dec[ix], 10) 

761 

762 def testSphericalFromCartesian(self): 

763 """ 

764 Note that xyz[i][j] is the ith component of the jth vector 

765 

766 Each column of xyz is a vector 

767 """ 

768 nsamples = 10 

769 radius = self.rng.random_sample(nsamples) * 10.0 

770 theta = self.rng.random_sample(nsamples) * np.pi - 0.5 * np.pi 

771 phi = self.rng.random_sample(nsamples) * 2.0 * np.pi 

772 

773 points = [] 

774 for ix in range(nsamples): 

775 vv = [radius[ix] * np.cos(theta[ix]) * np.cos(phi[ix]), 

776 radius[ix] * np.cos(theta[ix]) * np.sin(phi[ix]), 

777 radius[ix] * np.sin(theta[ix])] 

778 

779 points.append(vv) 

780 

781 points = np.array(points) 

782 lon, lat = utils.sphericalFromCartesian(points) 

783 for ix in range(nsamples): 

784 self.assertAlmostEqual(np.cos(lon[ix]), np.cos(phi[ix]), 5) 

785 self.assertAlmostEqual(np.sin(lon[ix]), np.sin(phi[ix]), 5) 

786 self.assertAlmostEqual(np.cos(lat[ix]), np.cos(theta[ix]), 5) 

787 self.assertAlmostEqual(np.sin(lat[ix]), np.sin(theta[ix]), 5) 

788 

789 # test passing in the points one at a time 

790 for pp, th, ph in zip(points, theta, phi): 

791 lon, lat = utils.sphericalFromCartesian(pp) 

792 self.assertAlmostEqual(np.cos(lon), np.cos(ph), 5) 

793 self.assertAlmostEqual(np.sin(lon), np.sin(ph), 5) 

794 self.assertAlmostEqual(np.cos(lat), np.cos(th), 5) 

795 self.assertAlmostEqual(np.sin(lat), np.sin(th), 5) 

796 

797 # test ra_dec_from_xyz <-> sphericalFromCartesian 

798 np.testing.assert_array_equal(utils.sphericalFromCartesian(points), 

799 utils._ra_dec_from_xyz(points[:, 0], points[:, 1], points[:, 2])) 

800 

801 # now, test passing one at a time 

802 for pp in points: 

803 np.testing.assert_array_equal(utils.sphericalFromCartesian(pp), 

804 utils._ra_dec_from_xyz(pp[0], pp[1], pp[2])) 

805 

806 def testCartesianFromSpherical(self): 

807 nsamples = 10 

808 theta = self.rng.random_sample(nsamples) * np.pi - 0.5 * np.pi 

809 phi = self.rng.random_sample(nsamples) * 2.0 * np.pi 

810 

811 points = [] 

812 for ix in range(nsamples): 

813 vv = [np.cos(theta[ix]) * np.cos(phi[ix]), 

814 np.cos(theta[ix]) * np.sin(phi[ix]), 

815 np.sin(theta[ix])] 

816 

817 points.append(vv) 

818 

819 points = np.array(points) 

820 lon, lat = utils.sphericalFromCartesian(points) 

821 outPoints = utils.cartesianFromSpherical(lon, lat) 

822 

823 for pp, oo in zip(points, outPoints): 

824 np.testing.assert_array_almost_equal(pp, oo, decimal=6) 

825 

826 # test passing in arguments as floats 

827 for ix, (ll, bb) in enumerate(zip(lon, lat)): 

828 xyz = utils.cartesianFromSpherical(ll, bb) 

829 self.assertIsInstance(xyz[0], np.float) 

830 self.assertIsInstance(xyz[1], np.float) 

831 self.assertIsInstance(xyz[2], np.float) 

832 self.assertAlmostEqual(xyz[0], outPoints[ix][0], 12) 

833 self.assertAlmostEqual(xyz[1], outPoints[ix][1], 12) 

834 self.assertAlmostEqual(xyz[2], outPoints[ix][2], 12) 

835 

836 # test _xyz_from_ra_dec <-> testCartesianFromSpherical 

837 np.testing.assert_array_equal(utils.cartesianFromSpherical(lon, lat), 

838 utils._xyz_from_ra_dec(lon, lat).transpose()) 

839 

840 def testHaversine(self): 

841 arg1 = 7.853981633974482790e-01 

842 arg2 = 3.769911184307751517e-01 

843 arg3 = 5.026548245743668986e+00 

844 arg4 = -6.283185307179586232e-01 

845 

846 output = utils.haversine(arg1, arg2, arg3, arg4) 

847 

848 self.assertAlmostEqual(output, 2.162615946398791955e+00, 10) 

849 

850 def testRotationMatrixFromVectors(self): 

851 v1 = np.zeros((3), dtype=float) 

852 v2 = np.zeros((3), dtype=float) 

853 v3 = np.zeros((3), dtype=float) 

854 

855 v1[0] = -3.044619987218469825e-01 

856 v2[0] = 5.982190522311925385e-01 

857 v1[1] = -5.473550908956383854e-01 

858 v2[1] = -5.573565912346714057e-01 

859 v1[2] = 7.795545496018386755e-01 

860 v2[2] = -5.757495946632366079e-01 

861 

862 output = utils.rotationMatrixFromVectors(v1, v2) 

863 

864 for i in range(3): 

865 for j in range(3): 

866 v3[i] += output[i][j]*v1[j] 

867 

868 for i in range(3): 

869 self.assertAlmostEqual(v3[i], v2[i], 7) 

870 

871 v1 = np.array([1.0, 1.0, 1.0]) 

872 self.assertRaises(RuntimeError, utils.rotationMatrixFromVectors, v1, v2) 

873 self.assertRaises(RuntimeError, utils.rotationMatrixFromVectors, v2, v1) 

874 

875class RotationTestCase(unittest.TestCase): 

876 

877 def setUp(self): 

878 self.sqrt2o2 = np.sqrt(2.0)/2.0 

879 self.x_vec = np.array([1.0, 0.0, 0.0]) 

880 self.y_vec = np.array([0.0, 1.0, 0.0]) 

881 self.z_vec = np.array([0.0, 0.0, 1.0]) 

882 self.px_py_vec = np.array([self.sqrt2o2, self.sqrt2o2, 0.0]) 

883 self.px_ny_vec = np.array([self.sqrt2o2, -self.sqrt2o2, 0.0]) 

884 self.nx_py_vec = np.array([-self.sqrt2o2, self.sqrt2o2, 0.0]) 

885 self.nx_ny_vec = np.array([-self.sqrt2o2, -self.sqrt2o2, 0.0]) 

886 

887 self.px_pz_vec = np.array([self.sqrt2o2, 0.0, self.sqrt2o2]) 

888 self.px_nz_vec = np.array([self.sqrt2o2, 0.0, -self.sqrt2o2]) 

889 self.nx_pz_vec = np.array([-self.sqrt2o2, 0.0, self.sqrt2o2]) 

890 self.nx_nz_vec = np.array([-self.sqrt2o2, 0.0, -self.sqrt2o2]) 

891 

892 self.py_pz_vec = np.array([0.0, self.sqrt2o2, self.sqrt2o2]) 

893 self.ny_pz_vec = np.array([0.0, -self.sqrt2o2, self.sqrt2o2]) 

894 self.py_nz_vec = np.array([0.0, self.sqrt2o2, -self.sqrt2o2]) 

895 self.ny_nz_vec = np.array([0.0, -self.sqrt2o2, -self.sqrt2o2]) 

896 

897 def test_rot_z(self): 

898 out = utils.rotAboutZ(self.x_vec, 0.5*np.pi) 

899 np.testing.assert_array_almost_equal(out, self.y_vec, decimal=10) 

900 out = utils.rotAboutZ(self.x_vec, -0.5*np.pi) 

901 np.testing.assert_array_almost_equal(out, -1.0*self.y_vec, decimal=10) 

902 out = utils.rotAboutZ(self.x_vec, np.pi) 

903 np.testing.assert_array_almost_equal(out, -1.0*self.x_vec, decimal=10) 

904 out = utils.rotAboutZ(self.x_vec, 0.25*np.pi) 

905 np.testing.assert_array_almost_equal(out, self.px_py_vec, decimal=10) 

906 out = utils.rotAboutZ(self.x_vec, -0.25*np.pi) 

907 np.testing.assert_array_almost_equal(out, self.px_ny_vec, decimal=10) 

908 out = utils.rotAboutZ(self.x_vec, 0.75*np.pi) 

909 np.testing.assert_array_almost_equal(out, self.nx_py_vec, decimal=10) 

910 out = utils.rotAboutZ(self.y_vec, 0.5*np.pi) 

911 np.testing.assert_array_almost_equal(out, -1.0*self.x_vec, decimal=10) 

912 out = utils.rotAboutZ(self.px_py_vec, 0.5*np.pi) 

913 np.testing.assert_array_almost_equal(out, self.nx_py_vec, decimal=10) 

914 out = utils.rotAboutZ(self.px_py_vec, 0.75*np.pi) 

915 np.testing.assert_array_almost_equal(out, -1.0*self.x_vec, decimal=10) 

916 

917 out = utils.rotAboutZ(np.array([self.x_vec, self.y_vec, self.nx_py_vec, self.nx_ny_vec]), 

918 -0.25*np.pi) 

919 np.testing.assert_array_almost_equal(out[0], self.px_ny_vec, decimal=10) 

920 np.testing.assert_array_almost_equal(out[1], self.px_py_vec, decimal=10) 

921 np.testing.assert_array_almost_equal(out[2], self.y_vec, decimal=10) 

922 np.testing.assert_array_almost_equal(out[3], -1.0*self.x_vec, decimal=10) 

923 

924 def test_rot_y(self): 

925 out = utils.rotAboutY(self.x_vec, 0.5*np.pi) 

926 np.testing.assert_array_almost_equal(out, -1.0*self.z_vec, decimal=10) 

927 out = utils.rotAboutY(self.z_vec, 0.5*np.pi) 

928 np.testing.assert_array_almost_equal(out, self.x_vec, decimal=10) 

929 out = utils.rotAboutY(self.px_pz_vec, 0.75*np.pi) 

930 np.testing.assert_array_almost_equal(out, -1.0*self.z_vec, decimal=10) 

931 out = utils.rotAboutY(self.px_pz_vec, -0.5*np.pi) 

932 np.testing.assert_array_almost_equal(out, self.nx_pz_vec, decimal=10) 

933 

934 out = utils.rotAboutY(np.array([self.px_pz_vec, self.nx_pz_vec, self.z_vec, self.x_vec]), 

935 -0.75*np.pi) 

936 

937 np.testing.assert_array_almost_equal(out[0], -1.0*self.x_vec, decimal=10) 

938 np.testing.assert_array_almost_equal(out[1], -1.0*self.z_vec, decimal=10) 

939 np.testing.assert_array_almost_equal(out[2], self.nx_nz_vec, decimal=10) 

940 np.testing.assert_array_almost_equal(out[3], self.nx_pz_vec, decimal=10) 

941 

942 def test_rot_x(self): 

943 out = utils.rotAboutX(self.y_vec, 0.5*np.pi) 

944 np.testing.assert_array_almost_equal(out, self.z_vec, decimal=10) 

945 out = utils.rotAboutX(self.y_vec, 0.75*np.pi) 

946 np.testing.assert_array_almost_equal(out, self.ny_pz_vec, decimal=10) 

947 out = utils.rotAboutX(self.z_vec, 0.5*np.pi) 

948 np.testing.assert_array_almost_equal(out, -1.0*self.y_vec, decimal=10) 

949 out = utils.rotAboutX(self.z_vec, -0.25*np.pi) 

950 np.testing.assert_array_almost_equal(out, self.py_pz_vec, decimal=10) 

951 out = utils.rotAboutX(self.py_nz_vec, -0.5*np.pi) 

952 np.testing.assert_array_almost_equal(out, self.ny_nz_vec, decimal=10) 

953 out = utils.rotAboutX(self.ny_nz_vec, 0.25*np.pi) 

954 np.testing.assert_array_almost_equal(out, -1.0*self.z_vec, decimal=10) 

955 

956 out = utils.rotAboutX(np.array([self.z_vec, self.py_pz_vec, self.ny_pz_vec, self.y_vec]), 

957 0.25*np.pi) 

958 

959 np.testing.assert_array_almost_equal(out[0], self.ny_pz_vec, decimal=10) 

960 np.testing.assert_array_almost_equal(out[1], self.z_vec, decimal=10) 

961 np.testing.assert_array_almost_equal(out[2], -1.0*self.y_vec, decimal=10) 

962 np.testing.assert_array_almost_equal(out[3], self.py_pz_vec, decimal=10) 

963 

964 

965 

966class MemoryTestClass(lsst.utils.tests.MemoryTestCase): 

967 pass 

968 

969if __name__ == "__main__": 969 ↛ 970line 969 didn't jump to line 970, because the condition on line 969 was never true

970 lsst.utils.tests.init() 

971 unittest.main()