Kstars

starhopper.cpp
1/*
2 SPDX-FileCopyrightText: 2010 Akarsh Simha <akarshsimha@gmail.com>
3
4 SPDX-License-Identifier: GPL-2.0-or-later
5*/
6
7#include "starhopper.h"
8
9#include "kstarsdata.h"
10#include "ksutils.h"
11#include "starcomponent.h"
12#include "skyobjects/starobject.h"
13
14#include <kstars_debug.h>
15
27
28QList<const StarObject *> StarHopper::computePath_const(const SkyPoint &src, const SkyPoint &dest, float fov_,
29 float maglim_, QStringList *metadata)
30{
31 fov = fov_;
32 maglim = maglim_;
33 start = &src;
34 end = &dest;
35
36 came_from.clear();
37 result_path.clear();
38
39 // Implements the A* search algorithm
40
46
47 qCDebug(KSTARS) << "StarHopper is trying to compute a path from source: " << src.ra().toHMSString()
48 << src.dec().toDMSString() << " to destination: " << dest.ra().toHMSString() << dest.dec().toDMSString()
49 << "; a starhop of " << src.angularDistanceTo(&dest).Degrees() << " degrees!";
50
51 oSet.append(&src);
52 g_score[&src] = 0;
53 h_score[&src] = src.angularDistanceTo(&dest).Degrees() / fov;
54 f_score[&src] = h_score[&src];
55
56 while (!oSet.isEmpty())
57 {
58 qCDebug(KSTARS) << "Next step";
59 // Find the node with the lowest f_score value
60 SkyPoint const *curr_node = nullptr;
61 double lowfscore = 1.0e8;
62
63 foreach (const SkyPoint *sp, oSet)
64 {
65 if (f_score[sp] < lowfscore)
66 {
67 lowfscore = f_score[sp];
68 curr_node = sp;
69 }
70 }
71 if (curr_node == nullptr)
72 continue;
73
74 qCDebug(KSTARS) << "Lowest fscore (vertex distance-plus-cost score) is " << lowfscore
75 << " with coords: " << curr_node->ra().toHMSString() << curr_node->dec().toDMSString()
76 << ". Considering this node now.";
77 if (curr_node == &dest || (curr_node != &src && h_score[curr_node] < 0.5))
78 {
79 // We are at destination
80 reconstructPath(came_from[curr_node]);
81 qCDebug(KSTARS) << "We've arrived at the destination! Yay! Result path count: " << result_path.count();
82
83 // Just a test -- try to print out useful instructions to the debug console. Once we make star hopper unexperimental, we should move this to some sort of a display
84 qCDebug(KSTARS) << "Star Hopping Directions: ";
85 const SkyPoint *prevHop = start;
86
87 for (auto &hopStar : result_path)
88 {
89 QString direction;
91 double pa; // should be 0 to 2pi
92
93 dms angDist = prevHop->angularDistanceTo(hopStar, &pa);
94
95 dms dmsPA;
96 dmsPA.setRadians(pa);
97 direction = KSUtils::toDirectionString(dmsPA);
98
99 if (metadata)
100 {
101 if (!patternNames.contains(hopStar))
102 {
103 spectralChar += hopStar->spchar();
104 starHopDirections = i18n(" Slew %1 degrees %2 to find an %3 star of mag %4 ", QString::number(angDist.Degrees(), 'f', 2),
105 direction, spectralChar, QString::number(hopStar->mag()));
106 qCDebug(KSTARS) << starHopDirections;
107 }
108 else
109 {
110 starHopDirections = i18n(" Slew %1 degrees %2 to find a(n) %3", QString::number(angDist.Degrees(), 'f', 2), direction,
111 patternNames.value(hopStar));
112 qCDebug(KSTARS) << starHopDirections;
113 }
114 metadata->append(starHopDirections);
115 qCDebug(KSTARS) << "Appended starHopDirections to metadata";
116 }
118 }
119 qCDebug(KSTARS) << " The destination is within a field-of-view";
120
121 return result_path;
122 }
123
124 oSet.removeOne(curr_node);
125 cSet.append(curr_node);
126
127 // FIXME: Make sense. If current node ---> dest distance is
128 // larger than src --> dest distance by more than 20%, don't
129 // even bother considering it.
130
131 if (h_score[curr_node] > h_score[&src] * 1.2)
132 {
133 qCDebug(KSTARS) << "Node under consideration has larger distance to destination (h-score) than start node! "
134 "Ignoring it.";
135 continue;
136 }
137
138 // Get the list of stars that are neighbours of this node
139 QList<StarObject *> neighbors;
140
141 // FIXME: Actually, this should be done in
142 // HorizontalToEquatorial, but we do it here because SkyPoint
143 // needs a lot of fixing to handle unprecessed and precessed,
144 // equatorial and horizontal coordinates nicely
145 SkyPoint *CurrentNode = const_cast<SkyPoint *>(curr_node);
146 //CurrentNode->deprecess(KStarsData::Instance()->updateNum());
147 CurrentNode->catalogueCoord(KStarsData::Instance()->updateNum()->julianDay());
148 //qCDebug(KSTARS) << "Calling starsInAperture";
149 StarComponent::Instance()->starsInAperture(neighbors, *curr_node, fov, maglim);
150 qCDebug(KSTARS) << "Choosing next node from a set of " << neighbors.count();
151 // Look for the potential next node
153
154 for (auto &nhd_node : neighbors)
155 {
156 if (cSet.contains(nhd_node))
157 continue;
158
159 // Compute the tentative g_score
161 bool tentative_better;
162 if (!oSet.contains(nhd_node))
163 {
164 oSet.append(nhd_node);
165 tentative_better = true;
166 }
168 tentative_better = true;
169 else
170 tentative_better = false;
171
173 {
174 came_from[nhd_node] = curr_node;
176 h_score[nhd_node] = nhd_node->angularDistanceTo(&dest).Degrees() / fov;
178 }
179 }
180 }
181 qCDebug(KSTARS) << "REGRET! Returning empty list!";
182 return QList<StarObject const *>(); // Return an empty QList
183}
184
185void StarHopper::reconstructPath(SkyPoint const *curr_node)
186{
187 if (curr_node != start)
188 {
189 StarObject const *s = dynamic_cast<StarObject const *>(curr_node);
190 Q_ASSERT(s);
191 result_path.prepend(s);
192 reconstructPath(came_from[s]);
193 }
194}
195
196float StarHopper::cost(const SkyPoint *curr, const SkyPoint *next)
197{
198 // This is a very heuristic method, that tries to produce a cost
199 // for each hop.
200
201 float netcost = 0;
202
203 // Convert 'next' into a StarObject
204 if (next == start)
205 {
206 // If the next hop is back to square one, junk it
207 return 1e8;
208 }
209 bool isThisTheEnd = (next == end);
210
211 float magcost, speccost;
212 magcost = speccost = 0;
213
214 if (!isThisTheEnd)
215 {
216 // We ought to be dealing with a star
217 StarObject const *nextstar = dynamic_cast<StarObject const *>(next);
219
220 // Test 1: How bright is the star?
221 magcost =
222 nextstar->mag() - 7.0 +
223 5 * log10(
224 fov); // The brighter, the better. FIXME: 8.0 is now an arbitrary reference to the average faint star. Should actually depend on FOV, something like log( FOV ).
225
226 // Test 2: Is the star strikingly red / yellow coloured?
227 QString SpType = nextstar->sptype();
228 char spclass = SpType.at(0).toLatin1();
229 speccost = (spclass == 'G' || spclass == 'K' || spclass == 'M') ? -0.3 : 0;
230 /*
231 // Test 3: Is the star in the general direction of the object?
232 // We use the cosine rule to find the angle between the hop direction, and the direction to destination
233 // a = side joining curr to end
234 // b = side joining curr to next
235 // c = side joining next to end
236 // C = angle between curr-next and curr-end
237 double sina, sinb, cosa, cosb;
238 curr->angularDistanceTo(dest).SinCos( &sina, &cosa );
239 curr->angularDistanceTo(next).SinCos( &sinb, &cosb );
240 double cosc = cos(next->angularDistanceTo(end).radians());
241 double cosC = ( cosc - cosa * cosb ) / (sina * sinb);
242 float dircost;
243 if( cosC < 0 ) // Wrong direction!
244 dircost = 1e8; // Some humongous number;
245 else
246 dircost = sqrt( 1 - cosC * cosC ) / cosC; // tan( C )
247 */
248 }
249
250 // Test 4: How far is the hop?
251 double distcost =
252 (curr->angularDistanceTo(next).Degrees() /
253 fov); // 1 "magnitude" incremental cost for 1 FOV. Is this even required, or is it just equivalent to halving our distance unit? I think it is required since the hop is not necessarily in the direction of the object -- asimha
254
255 // Test 5: How effective is the hop? [Might not be required with A*]
256 // double distredcost = -((src->angularDistanceTo( dest ).Degrees() - next->angularDistanceTo( dest ).Degrees()) * 60 / fov)*3; // 3 "magnitudes" for 1 FOV closer
257
258 // Test 6: Is the destination an asterism? Are there bright stars clustered nearby?
260 StarComponent::Instance()->starsInAperture(localNeighbors, *next, fov / 10, maglim + 1.0);
261 double stardensitycost = 1 - localNeighbors.count(); // -1 "magnitude" for every neighbouring star
262
263// Test 7: Identify star patterns
264
265#define RIGHT_ANGLE_THRESHOLD 0.05
266#define EQUAL_EDGE_THRESHOLD 0.025
267
268 double patterncost = 0;
270 if (!isThisTheEnd)
271 {
272 // We ought to be dealing with a star
273 StarObject const *nextstar = dynamic_cast<StarObject const *>(next);
275
276 float factor = 1.0;
277 while (factor <= 10.0)
278 {
279 localNeighbors.clear();
280 StarComponent::Instance()->starsInAperture(
281 localNeighbors, *next, fov / factor,
282 nextstar->mag() + 1.0); // Use a larger aperture for pattern identification; max 1.0 mag difference
283 foreach (StarObject *star, localNeighbors)
284 {
285 if (star == nextstar)
286 localNeighbors.removeOne(star);
287 else if (fabs(star->mag() - nextstar->mag()) > 1.0)
288 localNeighbors.removeOne(star);
289 } // Now, we should have a pruned list
290 factor += 1.0;
291 if (localNeighbors.size() == 2)
292 break;
293 }
294 factor -= 1.0;
295 if (localNeighbors.size() == 2)
296 {
297 patternName = i18n("triangle (of similar magnitudes)"); // any three stars form a triangle!
298 // Try to find triangles. Note that we assume that the standard Euclidian metric works on a sphere for small angles, i.e. the celestial sphere is nearly flat over our FOV.
300 double dRA1 = nextstar->ra().radians() - star1->ra().radians();
301 double dDec1 = nextstar->dec().radians() - star1->dec().radians();
302 double dist1sqr = dRA1 * dRA1 + dDec1 * dDec1;
303
305 double dRA2 = nextstar->ra().radians() - star2->ra().radians();
306 double dDec2 = nextstar->dec().radians() - star2->dec().radians();
307 double dist2sqr = dRA2 * dRA2 + dDec2 * dDec2;
308
309 // Check for right-angled triangles (without loss of generality, right angle is at this vertex)
310 if (fabs((dRA1 * dRA2 - dDec1 * dDec2) / sqrt(dist1sqr * dist2sqr)) < RIGHT_ANGLE_THRESHOLD)
311 {
312 // We have a right angled triangle! Give -3 magnitudes!
313 patterncost += -3;
314 patternName = i18n("right-angled triangle");
315 }
316
317 // Check for isosceles triangles (without loss of generality, this is the vertex)
318 if (fabs((dist1sqr - dist2sqr) / (dist1sqr)) < EQUAL_EDGE_THRESHOLD)
319 {
320 patterncost += -1;
321 patternName = i18n("isosceles triangle");
322 if (fabs((dRA2 * dDec1 - dRA1 * dDec2) / sqrt(dist1sqr * dist2sqr)) < RIGHT_ANGLE_THRESHOLD)
323 {
324 patterncost += -1;
325 patternName = i18n("straight line of 3 stars");
326 }
327 // Check for equilateral triangles
328 double dist3 = star1->angularDistanceTo(star2).radians();
329 double dist3sqr = dist3 * dist3;
330 if (fabs((dist3sqr - dist1sqr) / dist1sqr) < EQUAL_EDGE_THRESHOLD)
331 {
332 patterncost += -1;
333 patternName = i18n("equilateral triangle");
334 }
335 }
336 }
337 // TODO: Identify squares.
338 if (!patternName.isEmpty())
339 {
340 patternName += i18n(" within %1% of FOV of the marked star", (int)(100.0 / factor));
341 patternNames.insert(nextstar, patternName);
342 }
343 }
344
346 if (netcost < 0)
347 netcost = 0.1; // FIXME: Heuristics aren't supposed to be entirely random. This one is.
348 qCDebug(KSTARS) << "Mag cost: " << magcost << "; Spec Cost: " << speccost << "; Dist Cost: " << distcost
349 << "; Density cost: " << stardensitycost << "; Pattern cost: " << patterncost << "; Net cost: " << netcost
350 << "; Pattern: " << patternName;
351 return netcost;
352}
float mag() const
Definition skyobject.h:206
The sky coordinates of a point in the sky.
Definition skypoint.h:45
const CachingDms & dec() const
Definition skypoint.h:269
const CachingDms & ra() const
Definition skypoint.h:263
dms angularDistanceTo(const SkyPoint *sp, double *const positionAngle=nullptr) const
Computes the angular distance between two SkyObjects.
Definition skypoint.cpp:899
static StarComponent * Instance()
QList< StarObject * > * computePath(const SkyPoint &src, const SkyPoint &dest, float fov__, float maglim__, QStringList *metadata_=nullptr)
Computes path for Star Hop.
This is a subclass of SkyObject.
Definition starobject.h:33
An angle, stored as degrees, but expressible in many ways.
Definition dms.h:38
const double & Degrees() const
Definition dms.h:141
QString i18n(const char *text, const TYPE &arg...)
QAction * next(const QObject *recvr, const char *slot, QObject *parent)
void clear()
bool contains(const Key &key) const const
iterator insert(const Key &key, const T &value)
T value(const Key &key) const const
void append(QList< T > &&value)
void clear()
qsizetype count() const const
void prepend(parameter_type value)
const QChar at(qsizetype position) const const
QString number(double n, char format, int precision)
This file is part of the KDE documentation.
Documentation copyright © 1996-2024 The KDE developers.
Generated on Fri May 17 2024 11:48:27 by doxygen 1.10.0 written by Dimitri van Heesch, © 1997-2006

KDE's Doxygen guidelines are available online.