source: subversion/applications/rendering/gosmore/density.c @ 18450

Revision 18450, 6.0 KB checked in by nic, 4 years ago (diff)

Fix bug #2415

  • Property svn:executable set to *
Line 
1// Run this program with :
2// gcc -Wall -lm density.c && ./a.exe <density4.txt >density.sh
3// sudo apt-get install netpbm
4// for n in *.pnm; do pnmtopng <$n >${n%.pnm}.png; done
5//
6// The input to it is a list of OO Calc expressions (sums)
7// Only the ranges are extracted. The actual calculations are ignored.
8// First it will output the osmosis command for generating the extracts
9// Then it will output a list of cells that are covered by none of the bboxes
10// (if any). It is recommended to add ranges until all cells are covered.
11// At the same time it will create and image map (density.html) and the
12// associated images.
13#include <stdio.h>
14#include <string.h>
15#include <math.h>
16
17#include "density.xbm"
18
19#define DIM 1024
20int col[2], row[2], c, se = 0, cov[DIM][DIM], i, j, m, mat[DIM][DIM];
21int hshrink, vshrink, bcnt[2] = { 0, 0 };
22char block[2][98765], *bptr[2] = { block[0], block[1] };
23
24#define N2(x) ((x) ? (x) + 'A' - 1 : '0')
25int main (void)
26  {
27    char fname[30];
28    FILE *pnm, *html, *html2, *gosm = fopen ("bboxes.c", "w");
29    memset (cov, 0, sizeof (cov));
30    memset (mat, 0, sizeof (mat));
31    html= fopen ("density.html", "w");
32    fprintf (html, "<html>\n\
33<head>\n\
34<script language=\"JavaScript\">\n\
35     last = \"\";\n\
36     function mouseover (itemID)\n\
37     {\n\
38         if (last != itemID) {\n\
39            last = itemID;\n\
40            document.getElementById (\"jrcmap\").src = itemID + '.png';\n\
41         }\n\
42     }\n\
43</script>\n\
44</head>\n\
45<body>\n\
46<br><br>\n\
47<p align=center><b>Select a rectangle.<br> (Hover the mouse)<br><br>\n\
48<map name=\"region-map\" id=\"region-map\" class=\"map-areas\">\n");
49    while ((c = getchar ()) != EOF) {
50      if (c == '(') col[0] = col[1] = row[0] = row[1] = se = 0;
51      if ('A' <= c && c <= 'Z') col[se] = col[se] * 26 + c - 'A' + 1;
52      if ('0' <= c && c <= '9') row[se] = row[se] * 10 + c - '0';
53      if (c == ':') se = 1;
54      if (c == ')' && se-- == 1) {
55//      printf ("%4d %4d - %4d %4d\n", col[0], row[0], col[1], row[1]);
56/*      sprintf (fname, "%c%c%c%04d%c%c%c%04d.poly", N2 (col[0]/26/26),
57                 N2 (col[0]/26%26, N2 (col[0] % 26 + 1), row[0],
58                 N2 (col[1]/26/26), N2 (col[1]/26%26), N2 (col[1] % 26 + 1),
59                 row[1]); */
60        sprintf (fname, "%04d%04d%04d%04d.pnm", --col[0], --row[0],
61                 col[1], row[1]);
62        fprintf (gosm, "{ %d, %d, %d, %d },\n", col[0] - 512, row[0] - 512,
63          col[1] - 512, row[1] - 512);
64        pnm = fopen (fname, "w");
65        fprintf (pnm, "P2\n%d %d\n255\n", DIM, DIM);
66        for (i = 0; i < DIM; i++) {
67          for (j = 0; j < DIM; j++) {
68#define P1Eq2(x) ((x) == 0 ? 2 : (x) > 0 ? 1 : 0)
69            fprintf (pnm, "%d\n", (map[i * DIM + j] < 99 ? 255 : 192) -
70              (P1Eq2 (i - row[0]) + P1Eq2 (row[1] - 1 - i) +
71               P1Eq2 (j - col[0]) + P1Eq2 (col[1] - 1 - j) >= 5 ? 192 : 0));
72          }
73        }
74        fclose (pnm);
75       
76        strcpy (fname + 16, ".html");
77        html2 = fopen (fname, "w");
78        fprintf (html2, "<html><body><br><br>\n\
79<p align=center><b><a href=\"http://www.openstreetmap.org/\
80?minlon=%.5lf&minlat=%.5lf&maxlon=%.5lf&maxlat=%.5lf&box=yes\">View %.16s\
81</a><br><br>\n<a href=\"%.16s.zip\">Download %.16s.zip</a></body></html>\n",
82                col[0] * 360.0 / DIM - 180,
83          (atan (exp ((1 - row[1] / 512.0) * M_PI)) - M_PI / 4) / M_PI * 360,
84                  col[1] * 360.0 / DIM - 180,
85          (atan (exp ((1 - row[0] / 512.0) * M_PI)) - M_PI / 4) / M_PI * 360,
86                fname, fname, fname);
87        fclose (html2);
88       
89        hshrink = col[1] - col[0] > 60 ? 15 : (col[1] - col[0]) / 4;
90        vshrink = row[1] - row[0] > 60 ? 15 : (row[1] - row[0]) / 4;
91        for (i = row[0] + vshrink; i <= row[1] - vshrink; i++) {
92          mat[i][col[0] + hshrink] = mat[i][col[1] - hshrink] = 192;
93        }
94        for (j = col[0] + hshrink; j <= col[1] - hshrink; j++) {
95          mat[row[0] + vshrink][j] = mat[row[1] - vshrink][j] = 192;
96        }
97        fprintf (html, "<area shape=\"rect\" coords=\"%d,%d,%d,%d\" href=\"%.16s.html\" \n\
98  OnMouseOver=\"mouseover('%.16s')\" OnMouseOut=\"mouseover('map')\" />\n",
99                 col[0] + hshrink, row[0] + vshrink,
100                 col[1] - hshrink, row[1] - vshrink, fname, fname);
101//      printf (" --bp file=%s --wx %.16s.osm.bz2 \\\n",fname, fname);
102        bptr[col[1] > 418 ? 0 : 1] += sprintf (bptr[col[1] > 418 ? 0 : 1],
103          " \\\n --bb  idTrackerType=\"BitSet\" left=%.5lf right=%.5lf top=%.5lf bottom=%.5lf --wx %.16s.osm.gz",
104                col[0] * 360.0 / DIM - 180, col[1] * 360.0 / DIM - 180,
105          (atan (exp ((1 - row[0] / 512.0) * M_PI)) - M_PI / 4) / M_PI * 360,
106          (atan (exp ((1 - row[1] / 512.0) * M_PI)) - M_PI / 4) / M_PI * 360,
107                fname);
108        bcnt[col[1] > 418 ? 0 : 1]++;
109/*      poly = fopen (fname, "w");
110        fprintf (poly, "%.16s\n1\n", fname);
111        for (i = 0; i < 4; i++) {
112          fprintf (poly, "  %10.5lf  %10.5lf\n", col[i * (3 - i) / 2] * 360.0
113                   / DIM - 180.0,  (atan (exp ((1 - row[i / 2] / 512.0)
114                                         * M_PI)) - M_PI / 4) / M_PI * 360);
115        }
116        fprintf (poly, "END\nEND\n");
117        fclose (poly); */
118        for (i = col[0]; i < col[1]; i++) {
119          for (j = row[0]; j < row[1]; j++) cov[i][j] = 1;
120        }
121      }
122    }
123    printf ("bzcat planet-latest.osm.bz2 | ./osmosis-0.31/bin/osmosis\
124 --read-xml enableDateParsing=no file=/dev/stdin --tee %d \\\n\
125  --bb idTrackerType=\"BitSet\" left=-180 right=-33.04688\
126 top=85.05113 bottom=-90 --wx left.osm.gz%s\n\n\
127gunzip <left.osm.gz | ./osmosis-0.31/bin/osmosis\
128 --read-xml enableDateParsing=no file=/dev/stdin --tee %d%s\n\n\
129for n in 0*.osm.gz\n\
130do gunzip <$n | ./gosmore rebuild\n\
131   mv gosmore.pak master.pak\n\
132   gunzip <$n | ./gosmore rebuild -85 -179 85 179\n\
133   a/zip nroets.openhost.dk/htdocs/bbox/${n:0:16}.zip gosmore.pak\n\
134done\n\n", bcnt[0] + 1, block[0], bcnt[1], block[1]);
135    for (i = 0; i < DIM; i++) {
136      for (j = 0; j < DIM; j++) {
137        if (!cov[i][j]) printf ("%c%c%c%d\n", (i)/26/26+'A'-1,
138                                (i)/26%26+'A'-1, (i)%26+'A' , j+1);
139      }
140    }
141    fprintf (html, "<img src=\"map.png\" id=\"jrcmap\" "
142                          "usemap=\"#region-map\" />\n</body>\n</html>\n");
143    fclose (html);
144    pnm = fopen ("map.pnm", "w");
145    fprintf (pnm, "P2\n%d %d\n255\n", DIM, DIM);
146    for (i = 0; i < DIM; i++) {
147      for (j = 0; j < DIM; j++) {
148        fprintf (pnm, "%d\n", (map[i * DIM + j] < 99 ? 255 : 192) - mat[i][j]);
149      }
150    }
151    fclose (pnm);
152    return 0;
153  }
Note: See TracBrowser for help on using the repository browser.