func ExampleRing() {
	// Example 45.a, p. 320
	earth, err := pp.LoadPlanet(pp.Earth, "")
	if err != nil {
		fmt.Println(err)
		return
	}
	saturn, err := pp.LoadPlanet(pp.Saturn, "")
	if err != nil {
		fmt.Println(err)
		return
	}
	B, Bʹ, ΔU, P, a, b := saturnring.Ring(2448972.50068, earth, saturn)
	fmt.Printf("B  = %.3f\n", B*180/math.Pi)
	fmt.Printf("Bʹ = %.3f\n", Bʹ*180/math.Pi)
	fmt.Printf("ΔU = %.3f\n", ΔU*180/math.Pi)
	fmt.Printf("P  = %.3f\n", P*180/math.Pi)
	fmt.Printf("a  = %.2d\n", base.NewFmtAngle(a))
	fmt.Printf("b  = %.2d\n", base.NewFmtAngle(b))
	// Output:
	// B  = 16.442
	// Bʹ = 14.679
	// ΔU = 4.198
	// P  = 6.741
	// a  = 35″.87
	// b  = 10″.15
}

func ExampleAberration() {
	// Example 23.a, p. 152
	α := base.NewRA(2, 46, 11.331).Rad()
	δ := base.NewAngle(false, 49, 20, 54.54).Rad()
	jd := julian.CalendarGregorianToJD(2028, 11, 13.19)
	Δα2, Δδ2 := apparent.Aberration(α, δ, jd)
	fmt.Printf("%.3s  %.3s\n", base.NewFmtAngle(Δα2), base.NewFmtAngle(Δδ2))
	// Output:
	// 30.045″  6.697″
}

func ExampleNutation() {
	// Example 23.a, p. 152
	α := base.NewRA(2, 46, 11.331).Rad()
	δ := base.NewAngle(false, 49, 20, 54.54).Rad()
	jd := julian.CalendarGregorianToJD(2028, 11, 13.19)
	Δα1, Δδ1 := apparent.Nutation(α, δ, jd)
	fmt.Printf("%.3s  %.3s\n", base.NewFmtAngle(Δα1), base.NewFmtAngle(Δδ1))
	// Output:
	// 15.843″  6.217″
}

func ExampleNewAngle() {
	// Example negative values, p. 9.
	a := base.NewAngle(true, 13, 47, 22)
	fmt.Println(base.NewFmtAngle(a.Rad()))
	a = base.NewAngle(true, 0, 32, 41)
	// use # flag to force output of all three components
	fmt.Printf("%#s\n", base.NewFmtAngle(a.Rad()))
	// Output:
	// -13°47′22″
	// -0°32′41″
}

func ExampleStellar() {
	// Exercise, p. 119.
	day1 := 7.
	day5 := 27.
	r2 := []float64{
		base.NewRA(15, 3, 51.937).Rad(),
		base.NewRA(15, 9, 57.327).Rad(),
		base.NewRA(15, 15, 37.898).Rad(),
		base.NewRA(15, 20, 50.632).Rad(),
		base.NewRA(15, 25, 32.695).Rad(),
	}
	d2 := []float64{
		base.NewAngle(true, 8, 57, 34.51).Rad(),
		base.NewAngle(true, 9, 9, 03.88).Rad(),
		base.NewAngle(true, 9, 17, 37.94).Rad(),
		base.NewAngle(true, 9, 23, 16.25).Rad(),
		base.NewAngle(true, 9, 26, 01.01).Rad(),
	}
	jd := julian.CalendarGregorianToJD(1996, 2, 17)
	dt := jd - base.J2000
	dy := dt / base.JulianYear
	dc := dy / 100
	fmt.Printf("%.2f years\n", dy)
	fmt.Printf("%.4f century\n", dc)

	pmr := -.649 // sec/cen
	pmd := -1.91 // sec/cen
	r1 := base.NewRA(15, 17, 0.421+pmr*dc).Rad()
	// Careful with quick and dirty way of applying correction to seconds
	// component before converting to radians.  The dec here is negative
	// so correction must be subtracted.  Alternative, less error-prone,
	// way would be to convert both to radians, then add.
	d1 := base.NewAngle(true, 9, 22, 58.54-pmd*dc).Rad()
	fmt.Printf("α′ = %.3d, δ′ = %.2d\n",
		base.NewFmtRA(r1), base.NewFmtAngle(d1))

	day, dd, err := conjunction.Stellar(day1, day5, r1, d1, r2, d2)
	if err != nil {
		fmt.Println(err)
		return
	}
	fmt.Println(base.NewFmtAngle(dd))
	dInt, dFrac := math.Modf(day)
	fmt.Printf("1996 February %d at %s TD\n", int(dInt),
		base.NewFmtTime(dFrac*24*3600))

	// Output:
	// -3.87 years
	// -0.0387 century
	// α′ = 15ʰ17ᵐ0ˢ.446, δ′ = -9°22′58″.47
	// 3′38″
	// 1996 February 18 at 6ʰ36ᵐ55ˢ TD
}

func ExampleEclipticAtHorizon() {
	ε := 23.44 * math.Pi / 180
	φ := 51 * math.Pi / 180
	θ := 75 * math.Pi / 180
	λ1, λ2, I := parallactic.EclipticAtHorizon(ε, φ, θ)
	fmt.Println(base.NewFmtAngle(λ1))
	fmt.Println(base.NewFmtAngle(λ2))
	fmt.Println(base.NewFmtAngle(I))
	// Output:
	// 169°21′30″
	// 349°21′30″
	// 61°53′14″
}

func TestDiurnalPathAtHorizon(t *testing.T) {
	φ := 40 * math.Pi / 180
	ε := 23.44 * math.Pi / 180
	J := parallactic.DiurnalPathAtHorizon(0, φ)
	Jexp := math.Pi/2 - φ
	if math.Abs((J-Jexp)/Jexp) > 1e-15 {
		t.Fatal("0 dec:", base.NewFmtAngle(J))
	}
	J = parallactic.DiurnalPathAtHorizon(ε, φ)
	Jexp = base.NewAngle(false, 45, 31, 0).Rad()
	if math.Abs((J-Jexp)/Jexp) > 1e-3 {
		t.Fatal("solstace:", base.NewFmtAngle(J))
	}
}

func ExampleElements_Position() {
	// Example 33.b, p. 232.
	earth, err := pp.LoadPlanet(pp.Earth, "")
	if err != nil {
		fmt.Println(err)
		return
	}
	k := &elliptic.Elements{
		TimeP: julian.CalendarGregorianToJD(1990, 10, 28.54502),
		Axis:  2.2091404,
		Ecc:   .8502196,
		Inc:   11.94524 * math.Pi / 180,
		Node:  334.75006 * math.Pi / 180,
		ArgP:  186.23352 * math.Pi / 180,
	}
	j := julian.CalendarGregorianToJD(1990, 10, 6)
	α, δ, ψ := k.Position(j, earth)
	fmt.Printf("α = %.1d\n", base.NewFmtRA(α))
	fmt.Printf("δ = %.0d\n", base.NewFmtAngle(δ))
	fmt.Printf("ψ = %.2f\n", ψ*180/math.Pi)
	// Output:
	// α = 10ʰ34ᵐ14ˢ.2
	// δ = 19°9′31″
	// ψ = 40.51
}

func ExampleProperMotion3D() {
	// Example 21.d, p. 141.
	eqFrom := &coord.Equatorial{
		RA:  base.NewRA(6, 45, 8.871).Rad(),
		Dec: base.NewAngle(true, 16, 42, 57.99).Rad(),
	}
	mra := base.NewHourAngle(false, 0, 0, -0.03847)
	mdec := base.NewAngle(false, 0, 0, -1.2053)
	r := 2.64           // given in correct unit
	mr := -7.6 / 977792 // magic conversion factor
	eqTo := &coord.Equatorial{}
	fmt.Printf("Δr = %.9f, Δα = %.10f, Δδ = %.10f\n", mr, mra, mdec)
	for _, epoch := range []float64{1000, 0, -1000, -2000, -10000} {
		precess.ProperMotion3D(eqFrom, eqTo, 2000, epoch, r, mr, mra, mdec)
		fmt.Printf("%8.1f  %0.2d  %-0.1d\n", epoch,
			base.NewFmtRA(eqTo.RA), base.NewFmtAngle(eqTo.Dec))
	}
	// Output:
	// Δr = -0.000007773, Δα = -0.0000027976, Δδ = -0.0000058435
	//   1000.0  6ʰ45ᵐ47ˢ.16  -16°22′56″.0
	//      0.0  6ʰ46ᵐ25ˢ.09  -16°03′00″.8
	//  -1000.0  6ʰ47ᵐ02ˢ.67  -15°43′12″.3
	//  -2000.0  6ʰ47ᵐ39ˢ.91  -15°23′30″.6
	// -10000.0  6ʰ52ᵐ25ˢ.72  -12°50′06″.7
}

// Exercise, p. 136.
func TestPosition(t *testing.T) {
	eqFrom := &coord.Equatorial{
		base.NewRA(2, 31, 48.704).Rad(),
		base.NewAngle(false, 89, 15, 50.72).Rad(),
	}
	eqTo := &coord.Equatorial{}
	mα := base.NewHourAngle(false, 0, 0, 0.19877)
	mδ := base.NewAngle(true, 0, 0, 0.0152)
	for _, tc := range []struct {
		α, δ string
		jde  float64
	}{
		{"1 22 33.90", "88 46 26.18", base.BesselianYearToJDE(1900)},
		{"3 48 16.43", "89 27 15.38", base.JulianYearToJDE(2050)},
		{"5 53 29.17", "89 32 22.18", base.JulianYearToJDE(2100)},
	} {
		epochTo := base.JDEToJulianYear(tc.jde)
		precess.Position(eqFrom, eqTo, 2000.0, epochTo, mα, mδ)
		αStr := fmt.Sprintf("%.2x", base.NewFmtRA(eqTo.RA))
		δStr := fmt.Sprintf("%.2x", base.NewFmtAngle(eqTo.Dec))
		if αStr != tc.α {
			t.Fatal("got:", αStr, "want:", tc.α)
		}
		if δStr != tc.δ {
			t.Fatal(δStr)
		}
	}
}

func ExampleApparentLongitude() {
	// Example 25.a, p. 165.
	T := base.J2000Century(julian.CalendarGregorianToJD(1992, 10, 13))
	fmt.Println("λ:", base.NewFmtAngle(solar.ApparentLongitude(T)))
	// Output:
	// λ: 199°54′32″
}

func ExampleAngleError() {
	// Example p. 125.
	rδ := base.NewRA(5, 32, 0.40).Rad()
	dδ := base.NewAngle(true, 0, 17, 56.9).Rad()
	rε := base.NewRA(5, 36, 12.81).Rad()
	dε := base.NewAngle(true, 1, 12, 7.0).Rad()
	rζ := base.NewRA(5, 40, 45.52).Rad()
	dζ := base.NewAngle(true, 1, 56, 33.3).Rad()

	n, ω := line.AngleError(rδ, dδ, rε, dε, rζ, dζ)
	fmt.Printf("%.62s\n", base.NewFmtAngle(n))
	fmt.Println(base.NewFmtAngle(ω))
	// Output:
	// 7°31′
	// -5′24″
}

func TestPrecessor_Precess(t *testing.T) {
	// Exercise, p. 136.
	eqFrom := &coord.Equatorial{
		RA:  base.NewRA(2, 31, 48.704).Rad(),
		Dec: base.NewAngle(false, 89, 15, 50.72).Rad(),
	}
	mα := base.NewHourAngle(false, 0, 0, .19877)
	mδ := base.NewAngle(false, 0, 0, -.0152)
	epochs := []float64{
		base.JDEToJulianYear(base.B1900),
		2050,
		2100,
	}
	answer := []string{
		"α = 1ʰ22ᵐ33ˢ.90   δ = +88°46′26″.18",
		"α = 3ʰ48ᵐ16ˢ.43   δ = +89°27′15″.38",
		"α = 5ʰ53ᵐ29ˢ.17   δ = +89°32′22″.18",
	}
	eqTo := &coord.Equatorial{}
	for i, epochTo := range epochs {
		precess.Position(eqFrom, eqTo, 2000, epochTo, mα, mδ)
		if answer[i] != fmt.Sprintf("α = %0.2d   δ = %+0.2d",
			base.NewFmtRA(eqTo.RA), base.NewFmtAngle(eqTo.Dec)) {
			t.Fatal(i)
		}
	}
}

// p. 83
func TestLatDiff(t *testing.T) {
	φ0 := base.NewAngle(false, 45, 5, 46.36).Rad()
	diff := base.NewFmtAngle(globe.GeocentricLatitudeDifference(φ0))
	if f := fmt.Sprintf("%.2d", diff); f != "11′32″.73" {
		t.Fatal(f)
	}
}

func ExampleNutation() {
	// Example 22.a, p. 148.
	jd := julian.CalendarGregorianToJD(1987, 4, 10)
	Δψ, Δε := nutation.Nutation(jd)
	ε0 := nutation.MeanObliquity(jd)
	ε := ε0 + Δε
	fmt.Printf("%+.3d\n", base.NewFmtAngle(Δψ))
	fmt.Printf("%+.3d\n", base.NewFmtAngle(Δε))
	fmt.Printf("%.3d\n", base.NewFmtAngle(ε0))
	fmt.Printf("%.3d\n", base.NewFmtAngle(ε))
	// Output:
	// -3″.788
	// +9″.443
	// 23°26′27″.407
	// 23°26′36″.850
}

func ExampleApogeeParallax() {
	// Example 50.a, p. 357.
	p := apsis.ApogeeParallax(1988.75)
	fmt.Printf("%.3f\n", p*180/math.Pi*3600)
	fmt.Printf("%0.3d\n", base.NewFmtAngle(p))
	// Output:
	// 3240.679
	// 54′00″.679
}

// Second exercise, p. 110.
func TestMinSep(t *testing.T) {
	sep, err := angle.MinSep(jd1, jd3, r1, d1, r2, d2)
	if err != nil {
		t.Fatal(err)
	}
	answer := .5017 * math.Pi / 180 // on p. 111
	if math.Abs((sep-answer)/sep) > 1e-3 {
		t.Fatal(base.NewFmtAngle(sep))
	}
}

func ExampleApparentEquatorial() {
	// Example 25.a, p. 165.
	jde := julian.CalendarGregorianToJD(1992, 10, 13)
	α, δ := solar.ApparentEquatorial(jde)
	fmt.Printf("α: %.1d\n", base.NewFmtRA(α))
	fmt.Printf("δ: %d\n", base.NewFmtAngle(δ))
	// Output:
	// α: 13ʰ13ᵐ31ˢ.4
	// δ: -7°47′6″
}

// First exercise, p. 110.
func TestSep(t *testing.T) {
	r1 := base.NewRA(4, 35, 55.2).Rad()
	d1 := base.NewAngle(false, 16, 30, 33).Rad()
	r2 := base.NewRA(16, 29, 24).Rad()
	d2 := base.NewAngle(true, 26, 25, 55).Rad()
	d := angle.Sep(r1, d1, r2, d2)
	answer := base.NewAngle(false, 169, 58, 0).Rad()
	if math.Abs(d-answer) > 1e-4 {
		t.Fatal(base.NewFmtAngle(d))
	}
}

func ExampleSep() {
	// Example 17.a, p. 110.
	r1 := base.NewRA(14, 15, 39.7).Rad()
	d1 := base.NewAngle(false, 19, 10, 57).Rad()
	r2 := base.NewRA(13, 25, 11.6).Rad()
	d2 := base.NewAngle(true, 11, 9, 41).Rad()
	d := angle.Sep(r1, d1, r2, d2)
	fmt.Println(base.NewFmtAngle(d))
	// Output:
	// 32°47′35″
}