Odvození zdánlivého směru Slunce

Rodrigues’ rotation formula -- otáčí vektor v okolo vektoru k o uhel f.

Rot[v_,k_,f_]:=v Cos[f]+Cross[k,v]Sin[f]+k(k.v)(1-Cos[f]);

Zadám si standardní (pravotočivou) bázi R3, kde z je směr ke Slunci (k obloze), x +sever/-jih a y +západ/-východ

vec0={{1,0,0},{0,1,0},{0,0,1}};
k1={0,1,0};
k2={1,0,0};

A chci tu bázi transformovat na libovolný místo na Zemi v libovolnou roční dobu s libovolným sklonem rotace.
První rotace je kolem osy y o úhel θ -- nastavení zeměpisné šířky (+ odpovídá severní polokouli, - jižní)

vec1={Rot[vec0[[1]],k1,θ],Rot[vec0[[2]],k1,θ],Rot[vec0[[3]],k1,θ]}

Druhá rotace je kolem osy x o úhel ϕ -- nastavení zeměpisné délky, resp. denní doby (0° odpovídá poledni, + odpovídá rotaci východním směrem, - západním; tzn. zvětšující se úhel odpovídá rotace ze západu na východ -- tedy stejně jako se točí Země)

vec2=FullSimplify[{Rot[vec1[[1]],k2,ϕ],Rot[vec1[[2]],k2,ϕ],Rot[vec1[[3]],k2,ϕ]}]

Třetí je rotace kolem osy y o úhel α -- sklon osy rotace (0° odpovídá kolmosti k rovině oběhu, +úhel je příklon ke Slunci)

vec3=FullSimplify[{Rot[vec2[[1]],k1,-α],Rot[vec2[[2]],k1,-α],Rot[vec2[[3]],k1,-α]}]

A poslední je opět kolem osy x o úhel δ -- poloha na oběžné dráze (0° odpovídá “létu”)

vec4=FullSimplify[{Rot[vec3[[1]],k2,δ],Rot[vec3[[2]],k2,δ],Rot[vec3[[3]],k2,δ]}]

Nyní je třeba v nové bázi vyjádřit směr ke Slunci. Ve staré bázi to bylo (0,0,1).

nrule=FullSimplify[Solve[{n1 vec4[[1,1]]+n2 vec4[[2,1]]+n3 vec4[[3,1]]==0,n1 vec4[[1,2]]+n2 vec4[[2,2]]+n3 vec4[[3,2]]==0,n1 vec4[[1,3]]+n2 vec4[[2,3]]+n3 vec4[[3,3]]==1},{n1,n2,n3}]]

{n1,n2,n3}/.nrule/.ϕ->ϕ-δ

MatrixForm[Transpose[{n1,n2,n3}/.nrule]]

Vykreslím si na zkoušku trajektorii Slunce

g1=Show[
Graphics3D[Cylinder[{{0,0,0},{0,0,-0.01}},1]],
Graphics3D[{Opacity[0.1],Sphere[]}],
Graphics3D[Line[{{1,0,0},{-1,0,0}}]],
Graphics3D[Line[{{0,1,0},{0,-1,0}}]],Boxed->False
];
g2=ParametricPlot3D[{n1,n2,n3}/.{{n1->Cos[δ] (-Cos[θ] Sin[α]+Cos[α] Cos[ϕ] Sin[θ])-Sin[δ] Sin[θ] Sin[ϕ],n2->Cos[ϕ] Sin[δ]+Cos[α] Cos[δ] Sin[ϕ],n3->Cos[α] Cos[δ] Cos[θ] Cos[ϕ]+Cos[δ] Sin[α] Sin[θ]-Cos[θ] Sin[δ] Sin[ϕ]}}/.{α->π/8,δ->0,θ->π/4},{ϕ,0,2π}];
Show[g1,g2,PlotRange->{{-1,1},{-1,1},{0,1}}]

Zkouška poledníků a rovnoběžek

Kód stažený z internetu pro dělání kružnicových oblouků

ClearAll[splineCircle];
splineCircle[m_List,r_,angles_List: {0,2 π}]:=Module[{seg,ϕ,start,end,pts,w,k},{start,end}=Mod[angles//N,2 π];
If[end≤start,end+=2 π];
seg=Quotient[end-start//N,π/2];
ϕ=Mod[end-start//N,π/2];
If[seg==4,seg=3;ϕ=π/2];
pts=r RotationMatrix[start].#&/@Join[Take[{{1,0},{1,1},{0,1},{-1,1},{-1,0},{-1,-1},{0,-1}},2 seg+1],RotationMatrix[seg π/2].#&/@{{1,Tan[ϕ/2]},{Cos[ϕ],Sin[ϕ]}}];
If[Length[m]==2,pts=m+#&/@pts,pts=m+#&/@Transpose[Append[Transpose[pts],ConstantArray[0,Length[pts]]]]];
w=Join[Take[{1,1/Sqrt[2],1,1/Sqrt[2],1,1/Sqrt[2],1},2 seg+1],{Cos[ϕ/2],1}];
k=Join[{0,0,0},Riffle[#,#]&@Range[seg+1],{seg+1}];
BSplineCurve[pts,SplineDegree->2,SplineKnots->k,SplineWeights->w]]/;Length[m]==2||Length[m]==3

Ale mně stačí obyčejná kružnička:

splineCircle[{0,0,0},1]

Takže se inspiruju a utvořím jednodušší a rychlejší kód:

splineCircle2[center_,radius_]:=BSplineCurve[Table[center+{{1,0,0},{1,1,0},{0,1,0},{-1,1,0},{-1,0,0},{-1,-1,0},{0,-1,0},{1,-1,0},{1,0,0}}[[i]]radius,{i,1,9}],SplineDegree->2,SplineKnots->{0,0,0,1,1,2,2,3,3,4,4,4},SplineWeights->{1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1}]

A použiju na vykreslení poledníků a rovnoběžek

rovnobezky[num_]:=Table[splineCircle2[{0,0,Cos[angle]},Sin[angle]],{angle,(π/(num+1)),π,(π/(num+1))}];
poledniky[num_]:=Table[Rotate[Rotate[splineCircle2[{0,0,0},1],π/2,{1,0,0}],angle,{0,0,1}],{angle,(π/num),π,(π/num)}];
Graphics3D[Join[rovnobezky[5],poledniky[5]]]

Zkouška zeměkoulového grafu

Manipulate[
Show[
Graphics3D[
{
{LightCyan,Specularity[White,50],Glow[GrayLevel[0.2]],Sphere[{0,0,0}]},
Rotate[Rotate[
{
{Red,Cylinder[{{0,0,1.2},{0,0,-1.2}},0.05]},
Rotate[{Black,Join[rovnobezky[5],poledniky[3]]},ϕ-δ,{0,0,1}],
{PointSize->0.05,Red,Point[{Cos[ϕ-δ]Cos[θ],Sin[ϕ-δ]Cos[θ],Sin[θ]}]}
},
α,{0,1,0}],δ,{0,0,1}]
},
Lighting->{{"Directional",White,{{1,0,0},{0,0,0}}}},
ViewPoint->{0,0,2},
ViewVertical->{0,-1,0},
Boxed->False
],PlotRange->{{-1.2,1.2},{-1.2,1.2},{-1.2,1.2}}
],
{{α,0,Dynamic["Sklon osy rotace α = "<>ToString[Round[360α/(2π)]]<>"°"]},0,π/2},
{{δ,0,Dynamic["Poloha na oběžné dráze δ = "<>ToString[Round[360δ/(2π)]]<>"°"]},0,2π},
{{ϕ,0,Dynamic["Zeměpisná délka ϕ = "<>ToString[Round[360ϕ/(2π)]]<>"°"]},-π,π},
{{θ,0,Dynamic["Zeměpisná šířka θ = "<>ToString[Round[360θ/(2π)]]<>"°"]},-π/2,π/2},
Initialization:>{
n[pϕ_,pθ_,pα_,pδ_]:={Cos[pδ] (Cos[pθ] Sin[pα]-Cos[pα] Cos[pϕ] Sin[pθ])+Sin[pδ] Sin[pθ] Sin[pϕ],Cos[pϕ] Sin[pδ]+Cos[pα] Cos[pδ] Sin[pϕ],Cos[pα] Cos[pδ] Cos[pθ] Cos[pϕ]+Cos[pδ] Sin[pα] Sin[pθ]-Cos[pθ] Sin[pδ] Sin[pϕ]},
splineCircle2[center_,radius_]:=BSplineCurve[Table[center+{{1,0,0},{1,1,0},{0,1,0},{-1,1,0},{-1,0,0},{-1,-1,0},{0,-1,0},{1,-1,0},{1,0,0}}[[i]]radius,{i,1,9}],SplineDegree->2,SplineKnots->{0,0,0,1,1,2,2,3,3,4,4,4},SplineWeights->{1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1,(1/(Sqrt[2])),1}],
rovnobezky[num_]:=Table[splineCircle2[{0,0,Cos[angle]},Sin[angle]],{angle,(π/(num+1)),π,(π/(num+1))}],
poledniky[num_]:=Table[Rotate[Rotate[splineCircle2[{0,0,0},1],π/2,{1,0,0}],angle,{0,0,1}],{angle,(π/num),π,(π/num)}],
CardinalDirections={"S","SSV","SV","VSV","V","VJV","JV","JJV","J","JJZ","JZ","ZJZ","Z","ZSZ","SZ","SSZ","S"},
CardinalDirection[ang_]:=CardinalDirections[[Ceiling[(ang+π/16)/(π/8)]]]
}
]

Převedení úhlu na světové strany

CardinalDirections={"S","SSV","SV","VSV","V","VJV","JV","JJV","J","JJZ","JZ","ZJZ","Z","ZSZ","SZ","SSZ","S"}

CardinalDirections={"S","SSV","SV","VSV","V","VJV","JV","JJV","J","JJZ","JZ","ZJZ","Z","ZSZ","SZ","SSZ","S"}
CardinalDirection[ang_]:=CardinalDirections[[Ceiling[(ang+π/16)/(π/8)]]]

CardinalDirection[355*2π/360]

Kód finálního grafu

Manipulate[
Switch[presetRovnobezky,
1,{θ=0.;presetRovnobezky=0},
2,{θ=α;presetRovnobezky=0},
3,{θ=0.873;presetRovnobezky=0},
4,{θ=π/2.-α;presetRovnobezky=0},
5,{θ=π/2.;presetRovnobezky=0}];
Show[
Graphics3D[{
{Red,Cylinder[{{0,0,0},{0,0,0.3}},0.1]},
{Glow[GrayLevel[0.2]],Cylinder[{{0,0,0},{0,0,-0.01}},1]},
{Opacity[0.1],Sphere[]},
Line[{{0,0,0},n[ϕ-δ,θ,α,δ]}],
Line[{{1,0,0.01},{-1,0,0.01}}],
Line[{{0,1,0.01},{0,-1,0.01}}],
Text["S",{0.9,0,0.15}],
Text["J",{-0.9,0,0.15}],
Text["Z",{0,0.9,0.15}],
Text["V",{0,-0.9,0.15}],
{Yellow,Glow[Yellow],Sphere[n[ϕ-δ,θ,α,δ],0.1]}
},Lighting->{{"Directional",If[n[ϕ-δ,θ,α,δ][[3]]>0,White,Black],{n[ϕ-δ,θ,α,δ],{0,0,0}}}}
],
ParametricPlot3D[n[pϕ,θ,α,δ],{pϕ,0,2.π}],
Boxed->False,
PlotRange->{{-1.1,1.1},{-1.1,1.1},{0,1.05}},
ViewPoint->view1[[1]],
ViewVertical->view1[[2]],
ImageSize->{500,500},
Epilog->{
Inset[
Graphics3D[{
{Yellow,Glow[Yellow],Sphere[{0,0,0},0.2]},
{LightBlue,Sphere[{Cos[δ],Sin[δ],0},0.1]},
{Red,Rotate[Cylinder[{{Cos[δ],Sin[δ],-0.15},{Cos[δ],Sin[δ],0.15}},0.01],α,{0,-1,0},{Cos[δ],Sin[δ],0}]},
splineCircle2[{0,0,0},1],
Switch[view2[[3]],
1,Arrow[{{Cos[δ],Sin[δ],0.5},{Cos[δ],Sin[δ],0.15}}],
2,Arrow[{{Cos[δ]+0.5Sin[δ],Sin[δ]-0.5Cos[δ],0},{Cos[δ]+0.15Sin[δ],Sin[δ]-0.15Cos[δ],0}}],
3,Arrow[{{0.5Cos[δ],0.5Sin[δ],0},{0.85Cos[δ],0.85Sin[δ],0}}]
]
},
Boxed->False,
PlotRange->{{-1.2,1.2},{-1.2,1.2},{-0.2,0.5}},
ViewPoint->{0,-100,30},
Lighting->{{"Point",White,{0,0,0}}},
ImageSize->{200,200}
],
{1-view1[[3]],view1[[3]]}
],
Inset[
Show[
Graphics3D[
{
{LightCyan,Specularity[White,50],Glow[GrayLevel[0.2]],Sphere[{0,0,0}]},
Rotate[Rotate[
{
{Red,Cylinder[{{0,0,1.2},{0,0,-1.2}},0.05]},
Rotate[{Black,Join[rovnobezky[5],poledniky[3]]},ϕ-δ,{0,0,1}],
{Red,Glow[RGBColor[0.5,0,0]],Sphere[{Cos[ϕ-δ]Cos[θ],Sin[ϕ-δ]Cos[θ],Sin[θ]},0.1]}
},
α,{0,1,0}],δ,{0,0,1}]
},
Lighting->{{"Directional",White,{{1,0,0},{0,0,0}}}},
ViewPoint->view2[[1]],
ViewVertical->view2[[2]],
Boxed->False
],
PlotRange->{{-1.2,1.2},{-1.2,1.2},{-1.2,1.2}},
ImageSize->{150,150}
],
{view1[[3]],view1[[3]]}
]
}
],
Item[Style["Zdánlivý pohyb Slunce po obloze",Bold,16],Alignment->Center],
Delimiter,
Item[Style["Hlavní parametry\n",Bold,12],Alignment->Center],
Control[{{ϕ,0.0001,Dynamic["Čas "<>ToString[Floor[((ϕ+π)/(2π))24.]]<>"h "<>ToString[Round[60FractionalPart[((ϕ+π)/(2π))24.]]]<>"m"]},-π,π}],
Control[{{θ,0.873,Dynamic["Zeměpisná šířka θ = "<>ToString[Round[360θ/(2.π)]]<>"°"]},-π/2.,π/2.}],
Control[{{α,0.41,Dynamic["Sklon osy rotace α = "<>ToString[Round[360α/(2.π)]]<>"°"]},0,π/2.}],
Control[{{δ,0,Dynamic["Poloha na oběžné dráze δ = "<>ToString[Round[360δ/(2π)]]<>"°"]},0,2.π}],
Delimiter,
Item[Style["Nastavení pohledů\n",Bold,12],Alignment->Center],
Control[{{view1,{{-100,100,100},{0,0,1},0.8},"Pohled na horizontový graf:"},{
{{-100,100,100},{0,0,1},0.8}->"Iso",
{{0,0,∞},{1,0,0},0.6}->"Zvrchu",
{{-∞,0,0},{0,0,1},0.66}->"Z jihu",
{{0,∞,0},{0,0,1},0.66}->"Ze západu"}}],
Control[{{view2,{{0,∞,0},{0,-1,0},2},"Pohled na zeměkoulový graf:"},{
{{0,0,∞},{0,-1,0},1}->"Top",
{{0,∞,0},{0,-1,0},2}->"Left",
{{∞,0,0},{0,0,1},3}->"Front"}}],
Delimiter,
Item[Style["Presety jednotlivých parametrů\n",Bold,12],Alignment->Center],
Control[{{δ,0,"Roční doba: "},{
0->Tooltip["L slnv.","Letní slunovrat"],
π/2.->Tooltip["P rovnd.","Podzimní rovnodennost"],
1.*π->Tooltip["Z slnv.","Zimní slunovrat"],
3.π/2->Tooltip["J rovnd.","Jarní rovnodennost"]}}],
Control[{{presetRovnobezky,0,"Rovnoběžky: "},{
1->"Rovník",
2->"Obratník",
3->"Praha",
4->"Pol.kruh",
5->"Pól"}}],
Control[{{δ,π/8.,"Sklon osy: "},{
0->Tooltip["0°","Osa kolmá k rovině oběhu"],
0.41->Tooltip["Země","Sklon osy jako u Země"],
1.435->Tooltip["Uran","Sklon osy jako u Uranu"],
π/2.->Tooltip["90°","Osa ležící v rovině oběhu"]}}],
Delimiter,
Item[Style["Další informace\n",Bold,12],Alignment->Center],
Item[Dynamic["Výška Slunce nad obzorem je "<>ToString[Round[(360/(2π))ArcSin[n[ϕ-δ,θ,α,δ][[3]]],0.1]]<>"°"]],
Item[Dynamic["Azimut Slunce je "<>ToString[Round[180-(360/(2π))ArcTan[-n[ϕ-δ,θ,α,δ][[1]],-n[ϕ-δ,θ,α,δ][[2]]]]]<>"° ("<>CardinalDirection[π-ArcTan[-n[ϕ-δ,θ,α,δ][[1]],-n[ϕ-δ,θ,α,δ][[2]]]]<>")"]],
Item[Dynamic["Poledne je posunuto o "<>ToString[Round[24*60*(1/(2π))(ArcTan[-Cos[α]Cos[δ],+Sin[δ]]+δ-π),0.1]]<>" minut."]],
ControlPlacement->Left,
Initialization:>{
n[pϕ_,pθ_,pα_,pδ_]:={Cos[pδ] (Cos[pθ] Sin[pα]-Cos[pα] Cos[pϕ] Sin[pθ])+Sin[pδ] Sin[pθ] Sin[pϕ],Cos[pϕ] Sin[pδ]+Cos[pα] Cos[pδ] Sin[pϕ],Cos[pα] Cos[pδ] Cos[pθ] Cos[pϕ]+Cos[pδ] Sin[pα] Sin[pθ]-Cos[pθ] Sin[pδ] Sin[pϕ]},
splineCircle2[center_,radius_]:=BSplineCurve[Table[center+{{1.,0,0},{1.,1.,0},{0,1.,0},{-1.,1.,0},{-1.,0,0},{-1.,-1.,0},{0,-1.,0},{1.,-1.,0},{1.,0,0}}[[i]]radius,{i,1,9}],SplineDegree->2,SplineKnots->{0,0,0,1,1,2,2,3,3,4,4,4},SplineWeights->{1.,(1./(Sqrt[2])),1.,(1./(Sqrt[2])),1.,(1./(Sqrt[2])),1.,(1./(Sqrt[2])),1.}],
rovnobezky[num_]:=Table[splineCircle2[{0,0,Cos[angle]},Sin[angle]],{angle,1.*(π/(num+1)),1.*π,1.*(π/(num+1))}],
poledniky[num_]:=Table[Rotate[Rotate[splineCircle2[{0,0,0},1.],π/2.,{1,0,0}],angle,{0,0,1}],{angle,1.*(π/num),1.*π,1.*(π/num)}],
CardinalDirections={"S","SSV","SV","VSV","V","VJV","JV","JJV","J","JJZ","JZ","ZJZ","Z","ZSZ","SZ","SSZ","S"},
CardinalDirection[ang_]:=CardinalDirections[[Ceiling[(ang+π/16.)/(π/8.)]]]
}
]

Nějaké nedokončené výpočty...

FullSimplify[Solve[Cos[α]Cos[δ]Cos[θ]Cos[ϕ]+Cos[δ]Sin[α]Sin[θ]-Cos[θ]Sin[δ]Sin[ϕ]==0,ϕ]]

FullSimplify[Reduce[Cos[α]Cos[δ]Cos[θ]Cos[ϕ]+Cos[δ]Sin[α]Sin[θ]-Cos[θ]Sin[δ]Sin[ϕ]==0,ϕ]]

FullSimplify[{
ArcTan[-((Sec[θ] (Sec[θ] (Sqrt[((Cos[θ])^2) (((Cos[θ])^2)-((Cos[δ])^2) ((Sin[α])^2)) ((Sin[δ])^2)])+Cos[α] ((Cos[δ])^2) Sin[α] Sin[θ]))/(((Cos[α])^2) ((Cos[δ])^2)+((Sin[δ])^2))),((Cot[δ] (-Cos[α] ((Csc[δ])^2) ((Sec[θ])^2) (Sqrt[((Cos[θ])^2) (((Cos[θ])^2)-((Cos[δ])^2) ((Sin[α])^2)) ((Sin[δ])^2)])+Sin[α] Tan[θ]))/(1+((Cos[α])^2) ((Cot[δ])^2)))],
ArcTan[((((Sec[θ])^2) (Sqrt[((Cos[θ])^2) (((Cos[θ])^2)-((Cos[δ])^2) ((Sin[α])^2)) ((Sin[δ])^2)])-Cos[α] ((Cos[δ])^2) Sin[α] Tan[θ])/(((Cos[α])^2) ((Cos[δ])^2)+((Sin[δ])^2))),((Cot[δ] (Cos[α] ((Csc[δ])^2) ((Sec[θ])^2) (Sqrt[((Cos[θ])^2) (((Cos[θ])^2)-((Cos[δ])^2) ((Sin[α])^2)) ((Sin[δ])^2)])+Sin[α] Tan[θ]))/(1+((Cos[α])^2) ((Cot[δ])^2)))]
}]

NSolve[Cos[α]Cos[δ]Cos[θ]Cos[ϕ]+Cos[δ]Sin[α]Sin[θ]-Cos[θ]Sin[δ]Sin[ϕ]==0/.{α->23/360*2π,θ->50/360*2π,δ->0},ϕ]

angg=3;
ToString[Floor[(angg/(2π))24]]<>"h "<>ToString[Round[60FractionalPart[(angg/(2π))24]]]<>"m"

82./360*2π

vychod[pθ_,pα_,pδ_]:=ArcTan[-((Sec[pθ] (Sec[pθ] (Sqrt[((Cos[pθ])^2) (((Cos[pθ])^2)-((Cos[pδ])^2) ((Sin[pα])^2)) ((Sin[pδ])^2)])+Cos[pα] ((Cos[pδ])^2) Sin[pα] Sin[pθ]))/(((Cos[pα])^2) ((Cos[pδ])^2)+((Sin[pδ])^2))),((Cot[pδ] (-Cos[pα] ((Csc[pδ])^2) ((Sec[pθ])^2) (Sqrt[((Cos[pθ])^2) (((Cos[pθ])^2)-((Cos[pδ])^2) ((Sin[pα])^2)) ((Sin[pδ])^2)])+Sin[pα] Tan[pθ]))/(1+((Cos[pα])^2) ((Cot[pδ])^2)))],
zapad[pθ_,pα_,pδ_]:=ArcTan[((((Sec[pθ])^2) (Sqrt[((Cos[pθ])^2) (((Cos[pθ])^2)-((Cos[pδ])^2) ((Sin[pα])^2)) ((Sin[pδ])^2)])-Cos[pα] ((Cos[pδ])^2) Sin[pα] Tan[pθ])/(((Cos[pα])^2) ((Cos[pδ])^2)+((Sin[pδ])^2))),((Cot[pδ] (Cos[pα] ((Csc[pδ])^2) ((Sec[pθ])^2) (Sqrt[((Cos[pθ])^2) (((Cos[pθ])^2)-((Cos[pδ])^2) ((Sin[pα])^2)) ((Sin[pδ])^2)])+Sin[pα] Tan[pθ]))/(1+((Cos[pα])^2) ((Cot[pδ])^2)))]

FrontEndExecute[FrontEndToken["Save"]]