pd.if.org Git - isea/blob - transform.c

   1 #include <stdlib.h>
   2 #include <stdio.h>
   3 #include <math.h>
   4 #include <string.h>
   5
   6 #include "isea.h"
   7 #include "hex.h"
   8
   9 /*
  10  * options
  11  *
  12  * orientation
  13  * -l longitude
  14  * -z azimuth
  15  * -p latitude
  16  *
  17  * -I isea standard orientation
  18  * -D dymaxion orientation
  19  * -P vertex at pole
  20  *
  21  * -a aperture
  22  * -r resolution
  23  * -t h t d hex tri diamond
  24  *
  25  *  addressing
  26  *  -O q2dd seqnum interleave plane q2dd projtri vertex2dd
  27  *  -I
  28  *
  29  * operation
  30  * -g generate grid
  31  * -x transform
  32  * -b bin point values (average)
  33  * -e bin point presence
  34  * -G spacing resolution generate graticule
  35  */
  36
  37 void usage(void) {
  38         fprintf(stderr, "%s", "isea [-l <lon>] [-z <azimuth>] [-p <latitude>] [-I] [-D] [-P] [-a <aperture>] [-r <resolution>] [-t <topology>] [-A <addressing>] [-g] [-b] [-e] [-x]\n");
  39 }
  40
  41 void die(char *s) {
  42         fprintf(stderr, "%s\n", s);
  43         exit(EXIT_FAILURE);
  44 }
  45
  46 static int verbose = 0;
  47 static int precision = 10;
  48
  49 void options(int ac, char **av, struct isea_dgg *g) {
  50         int i = 1;
  51         char *optarg;
  52
  53         for (i=1;i<ac;i++) {
  54                 if (av[i][0] != '-') {
  55                         break;
  56                 }
  57                 switch (av[i][1]) {
  58                         case 'v':
  59                                 verbose = 1; break;
  60                         case 'I':
  61                                 isea_orient_isea(g); break;
  62                         case 'P':
  63                                 isea_orient_pole(g); break;
  64                         case 'D':
  65                                 isea_orient_dymax(g); break;
  66                         case 'R':
  67                                 g->radius = strtod(av[++i],NULL); break;
  68                         case 'l':
  69                                 g->o_lon = strtod(av[++i],NULL) * M_PI / 180.0; break;
  70                         case 'p':
  71                                 g->o_lat = strtod(av[++i],NULL) * M_PI / 180.0; break;
  72                         case 'z':
  73                                 g->o_az = strtod(av[++i],NULL) * M_PI / 180.0; break;
  74                         case 'd':
  75                                 precision = atoi(av[++i]); break;
  76                         case 'a':
  77                                 g->aperture = atoi(av[++i]); break;
  78                         case 'r':
  79                                 g->resolution = atoi(av[++i]); break;
  80                         case 't':
  81                                 switch (av[++i][0]) {
  82                                         case 'h':
  83                                                 g->topology = ISEA_HEXAGON;
  84                                                 break;
  85                                         case 't':
  86                                                 g->topology = ISEA_TRIANGLE;
  87                                                 break;
  88                                         case 'd':
  89                                                 g->topology = ISEA_DIAMOND;
  90                                                 break;
  91                                         default:
  92                                                 usage();
  93                                                 die("Invalid topology");
  94                                                 break;
  95                                 };
  96                                 break;
  97                         case 'O':
  98                                 if (av[i][2]) {
  99                                         optarg = &av[i][2];
 100                                 } else {
 101                                         optarg = av[++i];
 102                                 }
 103                                 if (!strcmp(optarg, "projtri")) {
 104                                         g->output = ISEA_PROJTRI;
 105                                 }
 106                                 else if (!strcmp(optarg, "vertex2dd")) {
 107                                         g->output = ISEA_VERTEX2DD;
 108                                 }
 109                                 else if (!strcmp(optarg, "q2dd")) {
 110                                         g->output = ISEA_Q2DD;
 111                                 }
 112                                 else if (!strcmp(optarg, "q2di")) {
 113                                         g->output = ISEA_Q2DI;
 114                                 }
 115                                 else if (!strcmp(optarg, "seqnum")) {
 116                                         g->output = ISEA_SEQNUM;
 117                                 }
 118                                 else if (!strcmp(optarg, "plane")) {
 119                                         g->output = ISEA_PLANE;
 120                                 } else {
 121                                         fprintf(stderr, "unknown output format: %s\n",optarg);
 122                                         die("exiting");
 123                                 }
 124                                 break;
 125                         default:
 126                                 usage();
 127                                 die("Unknown option");
 128                 }
 129         }
 130 }
 131
 132 int
 133 main(int ac, char *av[])
 134 {
 135         struct isea_dgg dgg;
 136         struct isea_geo in;
 137         struct isea_pt  out, coord;
 138         int             tri,calctri;
 139         int             graticule = 0;
 140         int     quad = 0;
 141
 142         char            line[1024];
 143         char           *rest, *nl;
 144
 145         isea_grid_init(&dgg);
 146         dgg.output = ISEA_PLANE;
 147         dgg.radius = ISEA_SCALE; /* 0.8301572857837594396028083; */
 148
 149         options(ac, av, &dgg);
 150
 151         if (verbose) {
 152                 fprintf(stderr, "isea lat = %lf, lon = %lf, az = %lf, output = %s\n",
 153                        dgg.o_lat * 180.0 / M_PI, dgg.o_lon * 180.0 / M_PI,
 154                        dgg.o_az * 180.0 / M_PI, dgg.output == ISEA_PLANE ? "plane" : "not plane");
 155         }
 156
 157         while (fgets(line, sizeof line, stdin)) {
 158                 line[sizeof line-1] = 0;
 159                 if (line[0] == '#') {
 160                         continue;
 161                 }
 162                 if (!strncmp(line, "ENDPOLY", 7)) {
 163                         printf("%s", line);
 164                         continue;
 165                 }
 166                 in.lon = strtod(line, &rest);
 167                 in.lat = strtod(rest, &rest);
 168
 169                 nl = strchr(rest, '\n');
 170                 if (nl) *nl = 0;
 171
 172                 in.lon *= M_PI / 180.0;
 173                 in.lat *= M_PI / 180.0;
 174
 175                 out = isea_forward(&dgg, &in);
 176                 tri = dgg.triangle;
 177                 quad = ((tri - 1) % 5) + ((tri - 1) / 10) * 5 + 1;
 178                 coord = out;
 179
 180                 if (dgg.output == ISEA_PLANE) {
 181                         calctri = isea_xy_triangle(&coord,
 182                                 dgg.output == ISEA_PLANE ? dgg.radius : 1.0);
 183 #if 0
 184                         if (tri != calctri) {
 185                                 fprintf(stderr,
 186                                         "got %d (%f, %f), expected %d\n",
 187                                         calctri, out.x, out.y, tri);
 188                         }
 189 #endif
 190                 }
 191
 192                 in.lon *= 180.0 / M_PI;
 193                 in.lat *= 180.0 / M_PI;
 194
 195                 switch (dgg.output) {
 196                 case ISEA_PLANE:
 197                         printf("%.*f %.*f%s\n", precision, out.x, precision, out.y, rest);
 198                         break;
 199                 case ISEA_PROJTRI:
 200                         printf("%d %.*f %.*f%s\n", tri-1, precision, out.x,
 201                                precision, out.y, rest);
 202                         break;
 203                 case ISEA_VERTEX2DD:
 204                         printf("%d %d %.*f %.*f%s\n", quad, tri-1, precision,
 205                                         out.x, precision, out.y, rest);
 206                         break;
 207                 case ISEA_Q2DD:
 208                         /* Same as above, we just don't print as much */
 209                         printf("%d %.*f %.*f%s\n", quad, precision,
 210                                         out.x, precision, out.y, rest);
 211                         break;
 212                 case ISEA_Q2DI:
 213                         printf("%d %.0lf %.0lf%s\n", dgg.quad, coord.x, coord.y, rest);
 214                         break;
 215                 case ISEA_SEQNUM:
 216                         quad = isea_ptdi(&dgg, tri, &out, &coord);
 217                         printf("%d%s\n", isea_disn(&dgg, quad, &coord), rest);
 218                         break;
 219                 }
 220                 fflush(stdout);
 221         }
 222
 223         if (graticule) {
 224
 225         }
 226
 227         exit(EXIT_SUCCESS);
 228 }