Python: Moved iau2000b and CalcMoon AddSol data into codegen phase.

Now the code generator knows how to optimize and inject these
two calculations into the Python code.
This commit is contained in:
Don Cross committed 2019-07-06 17:29:34 -04:00
1 parent 2920d983d3
commit d750c61337
5 files changed
+784 -3053

No files matched your search

+277
View File
@@ -31,6 +31,8 @@
#include "ephfile.h"
#define CG_MAX_LINE_LENGTH 200
#define MAX_DATA_PER_LINE 20
#define IAU_DATA_PER_ROW 11
static const double MJD_BASIS = 2400000.5;
@@ -645,6 +647,279 @@ fail:
return error;
}
static int ScanRealArray(
cg_context_t *context,
const char *filename,
int lnum,
char *line,
int numExpected,
double *data)
{
int i, t, len, inspace;
char *token[MAX_DATA_PER_LINE];
if (numExpected < 1 || numExpected > MAX_DATA_PER_LINE)
return LogError(context, "Invalid value for numExpected=%d\n", numExpected);
/* Split the line into space delimited tokens. */
len = strlen(line);
inspace = 1;
t = 0;
for (i=0; i < len; ++i)
{
if (line[i] == ' ' || line[i] == '\t' || line[i] == '\r' || line[i] == '\n')
{
if (!inspace)
{
line[i] = '\0'; /* terminate the previous token */
inspace = 1;
}
}
else
{
if (inspace)
{
/* we just found the beginning of a new token. */
if (t < IAU_DATA_PER_ROW)
{
token[t++] = &line[i];
inspace = 0;
}
else
return LogError(context, "ScanRealArray(%s %d): too many data on line.", filename, lnum);
}
}
}
/* Verify there are the correct number of tokens. */
if (t != numExpected)
return LogError(context, "ScanRealArray(%s %d): found %d data, but expected %d\n", filename, lnum, t, numExpected);
/* Parse each token as a floating point number. */
for (t=0; t < numExpected; ++t)
if (1 != sscanf(token[t], "%lf", &data[t]))
return LogError(context, "ScanRealArray(%s %d): invalid floating point token '%s'\n", filename, lnum, token[t]);
return 0; /* successful parse */
}
static int OptimizeConst(cg_context_t *context, char *buffer, size_t size, double c, const char *v)
{
int nprinted;
const char *op;
if (c == 0.0)
{
buffer[0] = '\0';
return 0;
}
if (c < 0.0)
{
op = " - ";
c *= -1.0;
}
else
op = " + ";
if (c == 1.0)
nprinted = snprintf(buffer, size, "%s%s", op, v);
else
nprinted = snprintf(buffer, size, "%s%0.1lf*%s", op, c, v);
if (nprinted >= (int)size)
return LogError(context, "OptimizeConst: print buffer overflowed.");
return 0;
}
static int OptimizeLinear(cg_context_t *context, char *buffer, size_t size, double a, double b)
{
int nprinted;
if (b == 0.0)
nprinted = snprintf(buffer, size, "%0.1lf", a);
else if (a == 0.0)
nprinted = snprintf(buffer, size, "%0.1lf*t", b);
else if (b < 0.0)
nprinted = snprintf(buffer, size, "%0.1lf - %0.1lf*t", a, -b);
else
nprinted = snprintf(buffer, size, "%0.1lf + %0.1lf*t", a, b);
if (nprinted >= size)
return LogError(context, "OptimizeLinear: print buffer overflowed.");
return 0;
}
static int OptIauPython(cg_context_t *context, const double *data)
{
static const char * const nv[] = {"el", "elp", "f", "d", "om"}; /* variable names */
int first;
double n;
int i;
const char *op;
char dotprod[200];
char term[40];
char linear[100];
char cpart[100];
int nprinted;
int lonevar;
fprintf(context->outfile, "\n"); /* must start on new line to maintain correct indentation in Python source. */
/* Optimize dot product of data[0]..data[4] with nv[]. */
first = 1;
dotprod[0] = '\0';
for (i=0; i<5; ++i)
{
n = data[i];
if (n != 0.0)
{
if (n < 0.0)
{
n *= -1.0;
op = first ? "-" : " - ";
}
else
op = first ? "" : " + ";
if (n == 1.0)
nprinted = snprintf(term, sizeof(term), "%s%s", op, nv[i]);
else
nprinted = snprintf(term, sizeof(term), "%s%0.1lf*%s", op, n, nv[i]);
if (nprinted >= sizeof(term))
return LogError(context, "Truncated iau2000b term.");
if (nprinted + strlen(dotprod) >= sizeof(dotprod))
return LogError(context, "Dot product overflow in iau2000b formula.");
strcat(dotprod, term);
first = 0;
}
}
/* Did we print a lone variable, e.g. "elp"? */
lonevar = 0;
for (i=0; i<5 && !lonevar; ++i)
if (!strcmp(dotprod, nv[i]))
lonevar = 1;
if (lonevar)
{
fprintf(context->outfile, " sarg = math.sin(%s)\n", dotprod);
fprintf(context->outfile, " carg = math.cos(%s)\n", dotprod);
}
else
{
fprintf(context->outfile, " arg = %s\n", dotprod);
fprintf(context->outfile, " sarg = math.sin(arg)\n");
fprintf(context->outfile, " carg = math.cos(arg)\n");
}
if (OptimizeLinear(context, linear, sizeof(linear), data[5], data[6])) return 1;
if (OptimizeConst(context, cpart, sizeof(cpart), data[7], "carg")) return 1;
fprintf(context->outfile, " dp += (%s)*sarg%s\n", linear, cpart);
if (OptimizeLinear(context, linear, sizeof(linear), data[8], data[9])) return 1;
if (OptimizeConst(context, cpart, sizeof(cpart), data[10], "sarg")) return 1;
fprintf(context->outfile, " de += (%s)*carg%s\n", linear, cpart);
fprintf(context->outfile, "\n");
return 0;
}
static int OptIauData(cg_context_t *context)
{
int error = 1;
int lnum;
FILE *infile;
const char *filename;
char line[100];
double data[IAU_DATA_PER_ROW];
filename = "model_data/iau2000b.txt";
infile = fopen(filename, "rt");
if (infile == NULL) goto fail;
lnum = 0;
while (fgets(line, sizeof(line), infile))
{
++lnum;
CHECK(ScanRealArray(context, filename, lnum, line, IAU_DATA_PER_ROW, data));
switch (context->language)
{
case CODEGEN_LANGUAGE_PYTHON:
CHECK(OptIauPython(context, data));
break;
default:
error = LogError(context, "OptIauData: Unsupported language %d", context->language);
goto fail;
}
}
error = 0;
fail:
if (infile == NULL)
error = LogError(context, "Cannot open input file: %s", filename);
else
fclose(infile);
return error;
}
static int OptAddSol(cg_context_t *context)
{
int nscanned;
double cl, cs, cg, cp, p, q, r, s;
const char *op;
nscanned = sscanf(context->args, "%lf , %lf , %lf , %lf , %lf , %lf , %lf , %lf", &cl, &cs, &cg, &cp, &p, &q, &r, &s);
if (nscanned != 8)
return LogError(context, "OptAddSol: invalid arguments: '%s'", context->args);
fprintf(context->outfile, "\n # AddSol(%s)\n", context->args);
op = "";
fprintf(context->outfile, " z = ");
if (p != 0.0)
{
fprintf(context->outfile, "ex[%0.0lf][1]", p);
op = " * ";
}
if (q != 0.0)
{
fprintf(context->outfile, "%sex[%0.0lf][2]", op, q);
op = " * ";
}
if (r != 0.0)
{
fprintf(context->outfile, "%sex[%0.0lf][3]", op, r);
op = " * ";
}
if (s != 0.0)
{
fprintf(context->outfile, "%sex[%0.0lf][4]", op, s);
}
fprintf(context->outfile, "\n");
if (cl != 0.0)
fprintf(context->outfile, " DLAM += %0.3lf * z.imag\n", cl);
if (cs != 0.0)
fprintf(context->outfile, " DS += %0.2lf * z.imag\n", cs);
if (cg != 0.0)
fprintf(context->outfile, " GAM1C += %0.3lf * z.real\n", cg);
if (cp != 0.0)
fprintf(context->outfile, " SINPI += %0.4lf * z.real\n", cp);
return 0;
}
static int LogError(const cg_context_t *context, const char *format, ...)
{
va_list v;
@@ -663,6 +938,8 @@ static const cg_directive_entry DirectiveTable[] =
{ "LIST_CHEBYSHEV", ListChebyshev },
{ "C_CHEBYSHEV", CChebyshev },
{ "DELTA_T", GenDeltaT },
{ "IAU_DATA", OptIauData },
{ "ADDSOL", OptAddSol },
{ NULL, NULL }
};
+77
View File
@@ -0,0 +1,77 @@
0 0 0 0 1 -172064161 -174666 33386 92052331 9086 15377
0 0 2 -2 2 -13170906 -1675 -13696 5730336 -3015 -4587
0 0 2 0 2 -2276413 -234 2796 978459 -485 1374
0 0 0 0 2 2074554 207 -698 -897492 470 -291
0 1 0 0 0 1475877 -3633 11817 73871 -184 -1924
0 1 2 -2 2 -516821 1226 -524 224386 -677 -174
1 0 0 0 0 711159 73 -872 -6750 0 358
0 0 2 0 1 -387298 -367 380 200728 18 318
1 0 2 0 2 -301461 -36 816 129025 -63 367
0 -1 2 -2 2 215829 -494 111 -95929 299 132
0 0 2 -2 1 128227 137 181 -68982 -9 39
-1 0 2 0 2 123457 11 19 -53311 32 -4
-1 0 0 2 0 156994 10 -168 -1235 0 82
1 0 0 0 1 63110 63 27 -33228 0 -9
-1 0 0 0 1 -57976 -63 -189 31429 0 -75
-1 0 2 2 2 -59641 -11 149 25543 -11 66
1 0 2 0 1 -51613 -42 129 26366 0 78
-2 0 2 0 1 45893 50 31 -24236 -10 20
0 0 0 2 0 63384 11 -150 -1220 0 29
0 0 2 2 2 -38571 -1 158 16452 -11 68
0 -2 2 -2 2 32481 0 0 -13870 0 0
-2 0 0 2 0 -47722 0 -18 477 0 -25
2 0 2 0 2 -31046 -1 131 13238 -11 59
1 0 2 -2 2 28593 0 -1 -12338 10 -3
-1 0 2 0 1 20441 21 10 -10758 0 -3
2 0 0 0 0 29243 0 -74 -609 0 13
0 0 2 0 0 25887 0 -66 -550 0 11
0 1 0 0 1 -14053 -25 79 8551 -2 -45
-1 0 0 2 1 15164 10 11 -8001 0 -1
0 2 2 -2 2 -15794 72 -16 6850 -42 -5
0 0 -2 2 0 21783 0 13 -167 0 13
1 0 0 -2 1 -12873 -10 -37 6953 0 -14
0 -1 0 0 1 -12654 11 63 6415 0 26
-1 0 2 2 1 -10204 0 25 5222 0 15
0 2 0 0 0 16707 -85 -10 168 -1 10
1 0 2 2 2 -7691 0 44 3268 0 19
-2 0 2 0 0 -11024 0 -14 104 0 2
0 1 2 0 2 7566 -21 -11 -3250 0 -5
0 0 2 2 1 -6637 -11 25 3353 0 14
0 -1 2 0 2 -7141 21 8 3070 0 4
0 0 0 2 1 -6302 -11 2 3272 0 4
1 0 2 -2 1 5800 10 2 -3045 0 -1
2 0 2 -2 2 6443 0 -7 -2768 0 -4
-2 0 0 2 1 -5774 -11 -15 3041 0 -5
2 0 2 0 1 -5350 0 21 2695 0 12
0 -1 2 -2 1 -4752 -11 -3 2719 0 -3
0 0 0 -2 1 -4940 -11 -21 2720 0 -9
-1 -1 0 2 0 7350 0 -8 -51 0 4
2 0 0 -2 1 4065 0 6 -2206 0 1
1 0 0 2 0 6579 0 -24 -199 0 2
0 1 2 -2 1 3579 0 5 -1900 0 1
1 -1 0 0 0 4725 0 -6 -41 0 3
-2 0 2 0 2 -3075 0 -2 1313 0 -1
3 0 2 0 2 -2904 0 15 1233 0 7
0 -1 0 2 0 4348 0 -10 -81 0 2
1 -1 2 0 2 -2878 0 8 1232 0 4
0 0 0 1 0 -4230 0 5 -20 0 -2
-1 -1 2 2 2 -2819 0 7 1207 0 3
-1 0 2 0 0 -4056 0 5 40 0 -2
0 -1 2 2 2 -2647 0 11 1129 0 5
-2 0 0 0 1 -2294 0 -10 1266 0 -4
1 1 2 0 2 2481 0 -7 -1062 0 -3
2 0 0 0 1 2179 0 -2 -1129 0 -2
-1 1 0 1 0 3276 0 1 -9 0 0
1 1 0 0 0 -3389 0 5 35 0 -2
1 0 2 0 0 3339 0 -13 -107 0 1
-1 0 2 -2 1 -1987 0 -6 1073 0 -2
1 0 0 0 2 -1981 0 0 854 0 0
-1 0 0 1 0 4026 0 -353 -553 0 -139
0 0 2 1 2 1660 0 -5 -710 0 -2
-1 0 2 4 2 -1521 0 9 647 0 4
-1 1 0 1 1 1314 0 0 -700 0 0
0 -2 2 -2 1 -1283 0 0 672 0 0
1 0 2 2 1 -1331 0 8 663 0 4
-2 0 2 2 2 1383 0 -2 -594 0 -2
-1 0 0 0 2 1405 0 4 -610 0 2
1 1 2 -2 2 1290 0 0 -556 0 0
File diff suppressed because it is too large. Load diff
File diff suppressed because it is too large. Load diff