<!--
    "Celestial Coordinates" calculates hour angle and right ascension
    Copyright (C) 2023  Matthias Liebendörfer

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <https://www.gnu.org/licenses/>.
-->

<phyphox version="1.12">
    <title>Celestial Coordinates</title>
    <category>Werkzeuge</category>
    <description>
        Launch the experiment and check in the tab "Location" whether
        date, time and location are correctly derived from the timer
        module and the GPS sensor.

        In the tab "RA-->HA" you may enter the right ascension (RA) of a
        celestial object and read out the hour angle (HA). The HA specifies
        the time after culmination of the object at your location.

        In the tab "POINTER" you get the declination and RA of the direction
        the long edge of your phone is pointing to. The orientation
        measurement is based on the magnetic field (and other sensors).
        Hence keep the phone reasonably far from magnetic materials.
        But even then the measurements may be inaccurate! Sometimes the
        accuracy improves with time or if you move the phone around and
        put it back into position.

        In this version of Celestial Coordinates the POINTER-functionality
        only works on android. It may work on ios as well using future
        versions of phyphox.
    </description>
    <data-containers>
        <container size="1">timeoffset</container>
        <container size="10">unixtime</container>
        <container size="1">jd</container>
        <container size="1">jdA</container>
        <container size="1">jdB</container>
        <container size="1">jdC</container>
        <container size="1">jdD</container>
        <container size="1">jdE</container>
        <container size="1">alpha</container>
        <container size="1">year</container>
        <container size="1">month</container>
        <container size="1">month1</container>
        <container size="1">month2</container>
        <container size="1">day</container>
        <container size="1">hour</container>
        <container size="1">minute</container>
        <container size="10">second</container>
        <container size="1">lat</container>
        <container size="1">lon</container>
        <container size="1">status</container>
        <container size="1">RAinh</container>
        <container size="1">RAinmin</container>
        <container size="1">tmp</container>
        <container size="1">GMST0</container>
        <container size="10">GMST1</container>
        <container size="10">LMST</container>
        <container size="10">HAdeg</container>
        <container size="10">HAh</container>
        <container size="10">HAmin</container>
        <container size="10">decl</container>
        <container size="10">az</container>
        <container size="10">RAdeg</container>
        <container size="10">RAh</container>
        <container size="10">RAmin</container>
        <container size="10">t</container>
        <container size="10">w</container>
        <container size="10">x</container>
        <container size="10">y</container>
        <container size="10">z</container>
        <container size="10">w1</container>
        <container size="10">x1</container>
        <container size="10">y1</container>
        <container size="10">z1</container>
        <container size="10">xvec</container>
        <container size="10">yvec</container>
        <container size="10">zvec</container>
        <container size="1" init="1">wq</container>
        <container size="1" init="0">xq</container>
        <container size="1" init="0">yq</container>
        <container size="1" init="0">zq</container>
    </data-containers>

    <!-- sensor input -->
    <input>
        <location>
            <output component="lat">lat</output>
            <output component="lon">lon</output>
            <output component="status">status</output>     
        </location>
        <sensor type="attitude" average="true" rate="1.0">
            <output component="x">x</output>
            <output component="y">y</output>
            <output component="z">z</output>
            <output component="abs">w</output>
        </sensor>
    </input>

    <!-- display -->
    <views>
        <view label="Location">
            <value label="Year" size="1" precision="0">
                <input>year</input>
            </value>
            <value label="Month" size="1" precision="0">
                <input>month</input>
            </value>
            <value label="Day" size="1" precision="0">
                <input>day</input>
            </value>
            <value label="UT Hour" unit="h" size="1" precision="0">
                <input>hour</input>
            </value>
            <value label="Minute" unit="min" size="1" precision="0">
                <input>minute</input>
            </value>
            <value label="Second" unit="s" size="1" precision="0">
                <input>second</input>
            </value>
            <value label="Latitude" unit="°" size="1" precision="3">
                <input>lat</input>
            </value>
            <value label="Longitude" unit="°" size="1" precision="3">
                <input>lon</input>
            </value>
            <value label="GPS Status" size="1" precision="0">
                <input>status</input>
            </value>
        </view>
        <view label="RA-->HA">
            <edit label="RA (hour)" signed="false" decimal="false" min="0" max="23">
                <output>RAinh</output>
            </edit>
            <edit label="RA (min)" signed="false" decimal="false" min="0" max="59">
                <output>RAinmin</output>
            </edit>
            <value color="red" label="HA (hour)" unit="h" size="3" precision="0">
                <input>HAh</input>
            </value>
            <value color="red" label="HA (min)" unit="'" size="3" precision="0">
                <input>HAmin</input>
            </value>
        </view>
        <view label="Pointer">
            <value color="red" label="Declination" unit="°" size="3" precision="1">
                <input>decl</input>
            </value>
            <value color="red" label="RA (hour)" unit="h" size="3" precision="0">
                <input>RAh</input>
            </value>
            <value color="red" label="RA (min)" unit="'" size="3" precision="0">
                <input>RAmin</input>
            </value>
        </view>
    </views>

    <!-- analysis of sensor data -->
    <analysis optimization="true">

        <!-- Calculate quaternion for rotation about x-axis -->
        <formula formula="cos(0.5*[1]/57.295779513)">
            <input clear="false">lat</input>
            <output clear="false">wq</output>
        </formula>
        <formula formula="-sin(0.5*[1]/57.295779513)">
            <input clear="false">lat</input>
            <output clear="false">xq</output>
        </formula>

        <!-- Multiply rotation about x-axis with attitude -->
        <formula formula="[1]*[5]-[2]*[6]-[3]*[7]-[4]*[8]">
            <input clear="false">wq</input>
            <input clear="false">xq</input>
            <input clear="false">yq</input>
            <input clear="false">zq</input>
            <input clear="false">w</input>
            <input clear="false">x</input>
            <input clear="false">y</input>
            <input clear="false">z</input>
            <output clear="false">w1</output>
        </formula>
        <formula formula="[1]*[6]+[2]*[5]+[3]*[8]-[4]*[7]">
            <input clear="false">wq</input>
            <input clear="false">xq</input>
            <input clear="false">yq</input>
            <input clear="false">zq</input>
            <input clear="false">w</input>
            <input clear="false">x</input>
            <input clear="false">y</input>
            <input clear="false">z</input>
            <output clear="false">x1</output>
        </formula>
        <formula formula="[1]*[7]+[3]*[5]+[4]*[6]-[2]*[8]">
            <input clear="false">wq</input>
            <input clear="false">xq</input>
            <input clear="false">yq</input>
            <input clear="false">zq</input>
            <input clear="false">w</input>
            <input clear="false">x</input>
            <input clear="false">y</input>
            <input clear="false">z</input>
            <output clear="false">y1</output>
        </formula>
        <formula formula="[1]*[8]+[4]*[5]+[2]*[7]-[3]*[6]">
            <input clear="false">wq</input>
            <input clear="false">xq</input>
            <input clear="false">yq</input>
            <input clear="false">zq</input>
            <input clear="false">w</input>
            <input clear="false">x</input>
            <input clear="false">y</input>
            <input clear="false">z</input>
            <output clear="false">z1</output>
        </formula>

        <!-- Apply combined quaternion to vector along local y-direction -->
        <formula formula="2*([2_]*[3_]-[1_]*[4_])">
            <input clear="false">w1</input>
            <input clear="false">x1</input>
            <input clear="false">y1</input>
            <input clear="false">z1</input>
            <output clear="false">xvec</output>
        </formula>
        <formula formula="[1_]*[1_]+[3_]*[3_]-[2_]*[2_]-[4_]*[4_]">
            <input clear="false">w1</input>
            <input clear="false">x1</input>
            <input clear="false">y1</input>
            <input clear="false">z1</input>
            <output clear="false">yvec</output>
        </formula>
        <formula formula="2*([1_]*[2_]+[3_]*[4_])">
            <input clear="false">w1</input>
            <input clear="false">x1</input>
            <input clear="false">y1</input>
            <input clear="false">z1</input>
            <output clear="false">zvec</output>
        </formula>

        <!-- Calculate declination of vector -->
        <formula formula="atan2([2_],sqrt([1_]*[1_]+[3_]*[3_]))*57.295779513">
            <input clear="false">xvec</input>
            <input clear="false">yvec</input>
            <input clear="false">zvec</input>
            <output clear="false">decl</output>
        </formula>

        <!-- Calculate azimuth of vector -->
        <formula formula="atan2(-[1_],[3_])*57.295779513">
            <input clear="false">xvec</input>
            <input clear="false">yvec</input>
            <input clear="false">zvec</input>
            <output clear="false">az</output>
        </formula>

        <!-- Calculate julian date from unix timestamp -->
        <timer linearTime="true">
            <output as="offset1970">timeoffset</output>
            <output as="out">t</output>
        </timer>
        <formula formula="[1]+[2_]">
            <input clear="false">timeoffset</input>
            <input clear="false">t</input>
            <output clear="false">unixtime</output>
        </formula>
        <!-- 1.1.1970 has julian date 2440587.5 -->
        <formula formula="2440587.5+floor([1_]/86400)">
            <input clear="false">unixtime</input>
            <output clear="false">jd</output>
        </formula>

        <!-- Calculate local time for display -->
        <!-- (https://de.wikipedia.org/wiki/Julianisches_Datum) -->
        <!-- jdZ = jd+0.5, jdF = 0 -->
        <formula formula="floor(([1]-1867215.75)/36524.25)">
            <input clear="false">jd</input>
            <output clear="false">alpha</output>
        </formula>
        <formula formula="[1]+1.5+[2]-floor([2]/4)">
            <input clear="false">jd</input>
            <input clear="false">alpha</input>
            <output clear="false">jdA</output>
        </formula>
        <formula formula="[1]+1524">
            <input clear="false">jdA</input>
            <output clear="false">jdB</output>
        </formula>
        <formula formula="floor(([1]-122.1)/365.25)">
            <input clear="false">jdB</input>
            <output clear="false">jdC</output>
        </formula>
        <formula formula="floor(365.25*[1])">
            <input clear="false">jdC</input>
            <output clear="false">jdD</output>
        </formula>
        <formula formula="floor(([1]-[2])/30.6001)">
            <input clear="false">jdB</input>
            <input clear="false">jdD</input>
            <output clear="false">jdE</output>
        </formula>
        <formula formula="[2]-[3]-floor(30.6001*[4])">
            <input clear="false">jd</input>
            <input clear="false">jdB</input>
            <input clear="false">jdD</input>
            <input clear="false">jdE</input>
            <output clear="false">day</output>
        </formula>
        <formula formula="[1]-1">
            <input clear="false">jdE</input>
            <output clear="false">month1</output>
        </formula>
        <formula formula="[1]-13">
            <input clear="false">jdE</input>
            <output clear="false">month2</output>
        </formula>
        <if less="true">
            <input clear="false">jdE</input>
            <input type="value">14</input>
            <input clear="false">month1</input>
            <input clear="false">month2</input>
            <output clear="false">month</output>
        </if>
        <formula formula="[1]-4716+floor(([2]-[3])/12)">
            <input clear="false">jdC</input>
            <input clear="false">jdE</input>
            <input clear="false">month</input>
            <output clear="false">year</output>
        </formula>
        <formula formula="floor([1_]%60)">
            <input clear="false">unixtime</input>
            <output clear="false">second</output>
        </formula>
        <formula formula="floor((([1_]-[2_])%3600)/60)">
            <input clear="false">unixtime</input>
            <input clear="false">second</input>
            <output clear="false">minute</output>
        </formula>
        <formula formula="floor((([1_]-[2_]-[3_]*60)%86400)/3600)">
            <input clear="false">unixtime</input>
            <input clear="false">second</input>
            <input clear="false">minute</input>
            <output clear="false">hour</output>
        </formula>

        <!-- Calculate local sideral time (https://de.wikipedia.org/wiki/Sternzeit) -->
        <formula formula="([1]-2451545.0)/36525">
            <input clear="false">jd</input>
            <output clear="false">tmp</output>
        </formula>
        <formula formula="24110.54841+[1]*(8640184.812866+[1]*(0.093104-[1]*0.0000062))">
            <input clear="false">tmp</input>
            <output clear="false">GMST0</output>
        </formula>
        <formula formula="[1]+([2_]%86400)*1.00273790935">
            <input clear="false">GMST0</input>
            <input clear="false">unixtime</input>
            <output clear="false">GMST1</output>
        </formula>
        <formula formula="[1_]/240+[2]">
            <input clear="false">GMST1</input>
            <input clear="false">lon</input>
            <output clear="false">LMST</output>
        </formula>

        <!-- Calculate hour angle (HA) -->
        <formula formula="([1_]-([2]*60+[3])/4)%360">
            <input clear="false">LMST</input>
            <input clear="false">RAinh</input>
            <input clear="false">RAinmin</input>
            <output clear="false">HAdeg</output>
        </formula>
        <formula formula="floor([1_]/15)">
            <input clear="false">HAdeg</input>
            <output clear="false">HAh</output>
        </formula>
        <formula formula="([1_]*4)%60">
            <input clear="false">HAdeg</input>
            <output clear="false">HAmin</output>
        </formula>

        <!-- Calculate right ascension (RA) -->
        <formula formula="([1_]-[2])%360">
            <input clear="false">LMST</input>
            <input clear="false">az</input>
            <output clear="false">RAdeg</output>
        </formula>
        <formula formula="floor([1_]/15)">
            <input clear="false">RAdeg</input>
            <output clear="false">RAh</output>
        </formula>
        <formula formula="([1_]*4)%60">
            <input clear="false">RAdeg</input>
            <output clear="false">RAmin</output>
        </formula>
    </analysis>

    <export>
        <set name="CelestialCoordinates">
            <data name="year">year</data>
            <data name="month">month</data>
            <data name="day">day</data>
            <data name="hour">hour</data>
            <data name="minute">minute</data>
            <data name="seconds">second</data>
            <data name="RAinh">RAinh</data>
            <data name="RAinmin">RAinmin</data>
            <data name="Hour Angle">HAdeg</data>
            <data name="Declination">decl</data>
            <data name="Right Ascension">RAdeg</data>
        </set>
    </export>
</phyphox>
