OpenCPN Partial API docs
Loading...
Searching...
No Matches
chartimg.cpp
1/**************************************************************************
2 * Copyright (C) 2015 by David S. Register *
3 * *
4 * This program is free software; you can redistribute it and/or modify *
5 * it under the terms of the GNU General Public License as published by *
6 * the Free Software Foundation; either version 2 of the License, or *
7 * (at your option) any later version. *
8 * *
9 * This program is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
12 * GNU General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU General Public License *
15 * along with this program; if not, see <https://www.gnu.org/licenses/>. *
16 **************************************************************************/
17
24#include <assert.h>
25
26// For compilers that support precompilation, includes "wx.h".
27#include <wx/wxprec.h>
28
29#ifndef WX_PRECOMP
30#include <wx/wx.h>
31#endif
32
33// Why are these not in wx/prec.h?
34#include <wx/dir.h>
35#include <wx/stream.h>
36#include <wx/wfstream.h>
37#include <wx/tokenzr.h>
38#include <wx/filename.h>
39#include <wx/image.h>
40#include <wx/fileconf.h>
41
43
44#include "config.h"
45#include "chartimg.h"
46#include "ocpn_pixel.h"
47
48#ifndef _MSC_VER
49#include <signal.h>
50#include <setjmp.h>
51#include "navutil.h"
52
53#define OCPN_USE_CONFIG 1
54#endif
55
56#ifndef _MSC_VER
57struct sigaction sa_all_chart;
58struct sigaction sa_all_previous;
59
60sigjmp_buf env_chart; // the context saved by sigsetjmp();
61
62void catch_signals_chart(int signo) {
63 switch (signo) {
64 case SIGSEGV:
65 siglongjmp(env_chart, 1); // jump back to the setjmp() point
66 break;
67
68 default:
69 break;
70 }
71}
72
73#endif
74
75// Missing from MSW include files
76#ifdef _MSC_VER
77typedef __int32 int32_t;
78typedef unsigned __int32 uint32_t;
79typedef __int64 int64_t;
80typedef unsigned __int64 uint64_t;
81#endif
82
83// ----------------------------------------------------------------------------
84// Random Prototypes
85// ----------------------------------------------------------------------------
86
87typedef struct {
88 float y;
89 float x;
90} MyFlPoint;
91
92bool G_FloatPtInPolygon(MyFlPoint *rgpts, int wnumpts, float x, float y);
93
94// ----------------------------------------------------------------------------
95// private classes
96// ----------------------------------------------------------------------------
97
98// ============================================================================
99// ThumbData implementation
100// ============================================================================
101
102ThumbData::ThumbData() { pDIBThumb = NULL; }
103
104ThumbData::~ThumbData() { delete pDIBThumb; }
105
106// ============================================================================
107// Palette implementation
108// ============================================================================
109opncpnPalette::opncpnPalette() {
110 // Index into palette is 1-based, so predefine the first entry as null
111 nFwd = 1;
112 nRev = 1;
113 FwdPalette = (int *)malloc(sizeof(int));
114 RevPalette = (int *)malloc(sizeof(int));
115 FwdPalette[0] = 0;
116 RevPalette[0] = 0;
117}
118
119opncpnPalette::~opncpnPalette() {
120 if (NULL != FwdPalette) free(FwdPalette);
121 if (NULL != RevPalette) free(RevPalette);
122}
123
124// ============================================================================
125// ChartBase implementation
126// ============================================================================
127ChartBase::ChartBase() {
128 m_depth_unit_id = DEPTH_UNIT_UNKNOWN;
129
130 pThumbData = new ThumbData;
131
132 m_global_color_scheme = GLOBAL_COLOR_SCHEME_RGB;
133
134 bReadyToRender = false;
135
136 Chart_Error_Factor = 0;
137
138 m_Chart_Scale = 10000; // a benign value
139 m_Chart_Skew = 0.0;
140
141 m_nCOVREntries = 0;
142 m_pCOVRTable = NULL;
143 m_pCOVRTablePoints = NULL;
144
145 m_nNoCOVREntries = 0;
146 m_pNoCOVRTable = NULL;
147 m_pNoCOVRTablePoints = NULL;
148
149 m_EdDate = wxInvalidDateTime;
150
151 m_lon_datum_adjust = 0.;
152 m_lat_datum_adjust = 0.;
153
154 m_projection = PROJECTION_MERCATOR; // default
155}
156
157ChartBase::~ChartBase() {
158 delete pThumbData;
159
160 // Free the COVR tables
161
162 for (unsigned int j = 0; j < (unsigned int)m_nCOVREntries; j++)
163 free(m_pCOVRTable[j]);
164
165 free(m_pCOVRTable);
166 free(m_pCOVRTablePoints);
167
168 // Free the No COVR tables
169
170 for (unsigned int j = 0; j < (unsigned int)m_nNoCOVREntries; j++)
171 free(m_pNoCOVRTable[j]);
172
173 free(m_pNoCOVRTable);
174 free(m_pNoCOVRTablePoints);
175}
176
177wxString ChartBase::GetHashKey() const {
178 wxString key = GetFullPath();
179 wxChar separator = wxFileName::GetPathSeparator();
180 for (unsigned int pos = 0; pos < key.size(); pos = key.find(separator, pos))
181 key.replace(pos, 1, "!");
182 return key;
183}
184
185/*
186int ChartBase::Continue_BackgroundHiDefRender(void)
187{
188 return BR_DONE_NOP; // signal "done, no refresh"
189}
190*/
191
192// ============================================================================
193// ChartDummy implementation
194// ============================================================================
195
196ChartDummy::ChartDummy() {
197 m_pBM = NULL;
198 m_ChartType = CHART_TYPE_DUMMY;
199 m_ChartFamily = CHART_FAMILY_UNKNOWN;
200 m_Chart_Scale = 22000000;
201
202 m_FullPath = "No Chart Available";
203 m_Description = m_FullPath;
204}
205
206ChartDummy::~ChartDummy() { delete m_pBM; }
207
208InitReturn ChartDummy::Init(const wxString &name, ChartInitFlag init_flags) {
209 return INIT_OK;
210}
211
212void ChartDummy::SetColorScheme(ColorScheme cs, bool bApplyImmediate) {}
213
214ThumbData *ChartDummy::GetThumbData(int tnx, int tny, float lat, float lon) {
215 return (ThumbData *)NULL;
216}
217
218bool ChartDummy::UpdateThumbData(double lat, double lon) { return FALSE; }
219
220bool ChartDummy::GetChartExtent(Extent *pext) {
221 pext->NLAT = 80;
222 pext->SLAT = -80;
223 pext->ELON = 179;
224 pext->WLON = -179;
225
226 return true;
227}
228
229bool ChartDummy::RenderRegionViewOnGL(const wxGLContext &glc,
230 const ViewPort &VPoint,
231 const OCPNRegion &RectRegion,
232 const LLRegion &Region) {
233 return true;
234}
235
236bool ChartDummy::RenderRegionViewOnDC(wxMemoryDC &dc, const ViewPort &VPoint,
237 const OCPNRegion &Region) {
238 return RenderViewOnDC(dc, VPoint);
239}
240
241bool ChartDummy::RenderViewOnDC(wxMemoryDC &dc, const ViewPort &VPoint) {
242 if (m_pBM && m_pBM->IsOk()) {
243 if ((m_pBM->GetWidth() != VPoint.pix_width) ||
244 (m_pBM->GetHeight() != VPoint.pix_height)) {
245 delete m_pBM;
246 m_pBM = NULL;
247 }
248 } else {
249 delete m_pBM;
250 m_pBM = NULL;
251 }
252
253 if (VPoint.pix_width && VPoint.pix_height) {
254 if (NULL == m_pBM)
255 m_pBM = new wxBitmap(VPoint.pix_width, VPoint.pix_height, -1);
256
257 dc.SelectObject(*m_pBM);
258
259 dc.SetBackground(*wxBLACK_BRUSH);
260 dc.Clear();
261 }
262
263 return true;
264}
265
266bool ChartDummy::AdjustVP(ViewPort &vp_last, ViewPort &vp_proposed) {
267 return false;
268}
269
270void ChartDummy::GetValidCanvasRegion(const ViewPort &VPoint,
271 OCPNRegion *pValidRegion) {
272 pValidRegion->Clear();
273 pValidRegion->Union(0, 0, 1, 1);
274}
275
276LLRegion ChartDummy::GetValidRegion() { return LLRegion(); }
277
278// ============================================================================
279// ChartGEO implementation
280// ============================================================================
281ChartGEO::ChartGEO() { m_ChartType = CHART_TYPE_GEO; }
282
283ChartGEO::~ChartGEO() {}
284
285InitReturn ChartGEO::Init(const wxString &name, ChartInitFlag init_flags) {
286#define BUF_LEN_MAX 4096
287
288 PreInit(name, init_flags, GLOBAL_COLOR_SCHEME_DAY);
289
290 char buffer[BUF_LEN_MAX];
291
292 ifs_hdr =
293 new wxFFileInputStream(name); // open the file as a read-only stream
294
295 m_filesize = wxFileName::GetSize(name);
296
297 if (!ifs_hdr->IsOk()) return INIT_FAIL_REMOVE;
298
299 int nPlypoint = 0;
300 Plypoint *pPlyTable = (Plypoint *)malloc(sizeof(Plypoint));
301
302 m_FullPath = name;
303 m_Description = m_FullPath;
304
305 wxFileName GEOFile(m_FullPath);
306
307 wxString Path;
308 Path = GEOFile.GetPath(wxPATH_GET_SEPARATOR | wxPATH_GET_VOLUME);
309
310 // Read the GEO file, extracting useful information
311
312 ifs_hdr->SeekI(0, wxFromStart); // rewind
313
314 Size_X = Size_Y = 0;
315
316 while ((ReadBSBHdrLine(ifs_hdr, &buffer[0], BUF_LEN_MAX)) != 0) {
317 wxString str_buf(buffer, wxConvUTF8);
318 if (!strncmp(buffer, "Bitmap", 6)) {
319 wxStringTokenizer tkz(str_buf, "=");
320 wxString token = tkz.GetNextToken();
321 if (token.IsSameAs("Bitmap", TRUE)) {
322 pBitmapFilePath = new wxString();
323
324 int i;
325 i = tkz.GetPosition();
326 pBitmapFilePath->Clear();
327 while (buffer[i]) {
328 pBitmapFilePath->Append(buffer[i]);
329 i++;
330 }
331 }
332 }
333
334 else if (!strncmp(buffer, "Scale", 5)) {
335 wxStringTokenizer tkz(str_buf, "=");
336 wxString token = tkz.GetNextToken();
337 if (token.IsSameAs("Scale", TRUE)) // extract Scale
338 {
339 int i;
340 i = tkz.GetPosition();
341 m_Chart_Scale = atoi(&buffer[i]);
342 }
343 }
344
345 else if (!strncmp(buffer, "Depth", 5)) {
346 wxStringTokenizer tkz(str_buf, "=");
347 wxString token = tkz.GetNextToken();
348 if (token.IsSameAs("Depth Units", FALSE)) // extract Depth Units
349 {
350 int i;
351 i = tkz.GetPosition();
352 wxString str(&buffer[i], wxConvUTF8);
353 m_DepthUnits = str.Trim();
354 }
355 }
356
357 else if (!strncmp(buffer, "Point", 5)) // Extract RefPoints
358 {
359 int i, xr, yr;
360 float ltr, lnr;
361 sscanf(&buffer[0], "Point%d=%f %f %d %d", &i, &lnr, &ltr, &yr, &xr);
362 pRefTable =
363 (Refpoint *)realloc(pRefTable, sizeof(Refpoint) * (nRefpoint + 1));
364 pRefTable[nRefpoint].xr = xr;
365 pRefTable[nRefpoint].yr = yr;
366 pRefTable[nRefpoint].latr = ltr;
367 pRefTable[nRefpoint].lonr = lnr;
368 pRefTable[nRefpoint].bXValid = 1;
369 pRefTable[nRefpoint].bYValid = 1;
370
371 nRefpoint++;
372
373 }
374
375 else if (!strncmp(buffer, "Vertex", 6)) {
376 int i;
377 float ltp, lnp;
378 sscanf(buffer, "Vertex%d=%f %f", &i, &ltp, &lnp);
379 Plypoint *tmp = pPlyTable;
380 pPlyTable =
381 (Plypoint *)realloc(pPlyTable, sizeof(Plypoint) * (nPlypoint + 1));
382 if (NULL == pPlyTable) {
383 free(tmp);
384 tmp = NULL;
385 } else {
386 pPlyTable[nPlypoint].ltp = ltp;
387 pPlyTable[nPlypoint].lnp = lnp;
388 nPlypoint++;
389 }
390 }
391
392 else if (!strncmp(buffer, "Date Pub", 8)) {
393 char date_string[40];
394 char date_buf[10];
395 sscanf(buffer, "Date Published=%s\r\n", &date_string[0]);
396 wxString date_wxstr(date_string, wxConvUTF8);
397 wxDateTime dt;
398 if (dt.ParseDate(date_wxstr)) // successful parse?
399 {
400 sprintf(date_buf, "%d", dt.GetYear());
401 } else {
402 sscanf(date_string, "%s", date_buf);
403 }
404 m_PubYear = wxString(date_buf, wxConvUTF8);
405 }
406
407 else if (!strncmp(buffer, "Skew", 4)) {
408 wxStringTokenizer tkz(str_buf, "=");
409 wxString token = tkz.GetNextToken();
410 if (token.IsSameAs("Skew Angle", FALSE)) // extract Skew Angle
411 {
412 int i;
413 i = tkz.GetPosition();
414 float fcs;
415 sscanf(&buffer[i], "%f,", &fcs);
416 m_Chart_Skew = fcs;
417 }
418 }
419
420 else if (!strncmp(buffer, "Latitude Offset", 15)) {
421 wxStringTokenizer tkz(str_buf, "=");
422 wxString token = tkz.GetNextToken();
423 if (token.IsSameAs("Latitude Offset", FALSE)) {
424 int i;
425 i = tkz.GetPosition();
426 float lto;
427 sscanf(&buffer[i], "%f,", &lto);
428 m_dtm_lat = lto;
429 }
430 }
431
432 else if (!strncmp(buffer, "Longitude Offset", 16)) {
433 wxStringTokenizer tkz(str_buf, "=");
434 wxString token = tkz.GetNextToken();
435 if (token.IsSameAs("Longitude Offset", FALSE)) {
436 int i;
437 i = tkz.GetPosition();
438 float lno;
439 sscanf(&buffer[i], "%f,", &lno);
440 m_dtm_lon = lno;
441 }
442 }
443
444 else if (!strncmp(buffer, "Datum", 5)) {
445 wxStringTokenizer tkz(str_buf, "=");
446 wxString token = tkz.GetNextToken();
447 if (token.IsSameAs("Datum", FALSE)) {
448 token = tkz.GetNextToken();
449 m_datum_str = token;
450 }
451 }
452
453 else if (!strncmp(buffer, "Name", 4)) {
454 wxStringTokenizer tkz(str_buf, "=");
455 wxString token = tkz.GetNextToken();
456 if (token.IsSameAs("Name", FALSE)) // Name
457 {
458 int i;
459 i = tkz.GetPosition();
460 m_Name.Clear();
461 while (isprint(buffer[i]) && (i < 80)) m_Name.Append(buffer[i++]);
462 }
463 }
464 } // while
465
466 // Extract the remaining data from .NOS Bitmap file
467 ifs_bitmap = NULL;
468
469 // Something wrong with the .geo file, there is no Bitmap reference
470 // This is where the arbitrarily bad file is caught, such as
471 // a file with.GEO extension that is not really a chart
472
473 if (pBitmapFilePath == NULL) {
474 free(pPlyTable);
475 return INIT_FAIL_REMOVE;
476 }
477
478 wxString NOS_Name(*pBitmapFilePath); // take a copy
479
480 wxDir target_dir(Path);
481 wxArrayString file_array;
482 int nfiles = wxDir::GetAllFiles(Path, &file_array);
483 int ifile;
484
485 pBitmapFilePath->Prepend(Path);
486
487 wxFileName NOS_filename(*pBitmapFilePath);
488 if (!NOS_filename.FileExists()) {
489 // File as fetched verbatim from the .geo file doesn't exist.
490 // Try all possible upper/lower cases
491 // Extract the filename and extension
492 wxString fname(NOS_filename.GetName());
493 wxString fext(NOS_filename.GetExt());
494
495 // Try all four combinations, the hard way
496 // case 1
497 fname.MakeLower();
498 fext.MakeLower();
499 NOS_filename.SetName(fname);
500 NOS_filename.SetExt(fext);
501
502 if (NOS_filename.FileExists()) goto found_uclc_file;
503
504 // case 2
505 fname.MakeLower();
506 fext.MakeUpper();
507 NOS_filename.SetName(fname);
508 NOS_filename.SetExt(fext);
509
510 if (NOS_filename.FileExists()) goto found_uclc_file;
511
512 // case 3
513 fname.MakeUpper();
514 fext.MakeLower();
515 NOS_filename.SetName(fname);
516 NOS_filename.SetExt(fext);
517
518 if (NOS_filename.FileExists()) goto found_uclc_file;
519
520 // case 4
521 fname.MakeUpper();
522 fext.MakeUpper();
523 NOS_filename.SetName(fname);
524 NOS_filename.SetExt(fext);
525
526 if (NOS_filename.FileExists()) goto found_uclc_file;
527
528 // Search harder
529
530 for (ifile = 0; ifile < nfiles; ifile++) {
531 wxString file_up = file_array[ifile];
532 file_up.MakeUpper();
533
534 wxString target_up = *pBitmapFilePath;
535 target_up.MakeUpper();
536
537 if (file_up.IsSameAs(target_up)) {
538 NOS_filename.Clear();
539 NOS_filename.Assign(file_array[ifile]);
540 goto found_uclc_file;
541 }
542 }
543
544 free(pPlyTable);
545 return INIT_FAIL_REMOVE; // not found at all
546
547 found_uclc_file:
548
549 delete pBitmapFilePath; // fix up the member element
550 pBitmapFilePath = new wxString(NOS_filename.GetFullPath());
551 }
552 ifss_bitmap =
553 new wxFFileInputStream(*pBitmapFilePath); // open the bitmap file
554 ifs_bitmap = new wxBufferedInputStream(*ifss_bitmap);
555
556 if (!ifss_bitmap->IsOk()) {
557 free(pPlyTable);
558 return INIT_FAIL_REMOVE;
559 }
560
561 while ((ReadBSBHdrLine(ifss_bitmap, &buffer[0], BUF_LEN_MAX)) != 0) {
562 wxString str_buf(buffer, wxConvUTF8);
563
564 if (!strncmp(buffer, "NOS", 3)) {
565 wxStringTokenizer tkz(str_buf, ",=");
566 while (tkz.HasMoreTokens()) {
567 wxString token = tkz.GetNextToken();
568 if (token.IsSameAs("RA", TRUE)) // extract RA=x,y
569 {
570 int i;
571 tkz.GetNextToken();
572 tkz.GetNextToken();
573 i = tkz.GetPosition();
574 Size_X = atoi(&buffer[i]);
575 wxString token = tkz.GetNextToken();
576 i = tkz.GetPosition();
577 Size_Y = atoi(&buffer[i]);
578 } else if (token.IsSameAs("DU", TRUE)) // extract DU=n
579 {
580 token = tkz.GetNextToken();
581 long temp_du;
582 if (token.ToLong(&temp_du)) m_Chart_DU = temp_du;
583 }
584 }
585
586 }
587
588 else if (!strncmp(buffer, "RGB", 3))
589 CreatePaletteEntry(buffer, COLOR_RGB_DEFAULT);
590
591 else if (!strncmp(buffer, "DAY", 3))
592 CreatePaletteEntry(buffer, DAY);
593
594 else if (!strncmp(buffer, "DSK", 3))
595 CreatePaletteEntry(buffer, DUSK);
596
597 else if (!strncmp(buffer, "NGT", 3))
598 CreatePaletteEntry(buffer, NIGHT);
599
600 else if (!strncmp(buffer, "NGR", 3))
601 CreatePaletteEntry(buffer, NIGHTRED);
602
603 else if (!strncmp(buffer, "GRY", 3))
604 CreatePaletteEntry(buffer, GRAY);
605
606 else if (!strncmp(buffer, "PRC", 3))
607 CreatePaletteEntry(buffer, PRC);
608
609 else if (!strncmp(buffer, "PRG", 3))
610 CreatePaletteEntry(buffer, PRG);
611 }
612
613 // Validate some of the header data
614 if (Size_X <= 0 || Size_Y <= 0) {
615 free(pPlyTable);
616 return INIT_FAIL_REMOVE;
617 }
618
619 if (nPlypoint < 3) {
620 wxString msg(" Chart File contains less than 3 PLY points: ");
621 msg.Append(m_FullPath);
622 wxLogMessage(msg);
623 free(pPlyTable);
624
625 return INIT_FAIL_REMOVE;
626 }
627
628 if (m_datum_str.IsEmpty()) {
629 wxString msg(" Chart datum not specified on chart ");
630 msg.Append(m_FullPath);
631 wxLogMessage(msg);
632 wxLogMessage(" Default datum (WGS84) substituted.");
633
634 // return INIT_FAIL_REMOVE;
635 } else {
636 char d_str[100];
637 strncpy(d_str, m_datum_str.mb_str(), 99);
638 d_str[99] = 0;
639
640 int datum_index = GetDatumIndex(d_str);
641
642 if (datum_index < 0) {
643 wxString msg(" Chart datum {");
644 msg += m_datum_str;
645 msg += "} invalid on chart ";
646 msg.Append(m_FullPath);
647 wxLogMessage(msg);
648 wxLogMessage(" Default datum (WGS84) substituted.");
649
650 datum_index = DATUM_INDEX_WGS84;
651 }
652 m_datum_index = datum_index;
653 }
654
655 // Convert captured plypoint information into chart COVR structures
656 m_nCOVREntries = 1;
657 m_pCOVRTablePoints = (int *)malloc(sizeof(int));
658 *m_pCOVRTablePoints = nPlypoint;
659 m_pCOVRTable = (float **)malloc(sizeof(float *));
660 *m_pCOVRTable = (float *)malloc(nPlypoint * 2 * sizeof(float));
661 memcpy(*m_pCOVRTable, pPlyTable, nPlypoint * 2 * sizeof(float));
662
663 free(pPlyTable);
664
665 if (!SetMinMax()) return INIT_FAIL_REMOVE; // have to bail here
666
667 AnalyzeSkew();
668
669 if (init_flags == HEADER_ONLY) return INIT_OK;
670
671 // Advance to the data
672 char c;
673 if ((c = ifs_bitmap->GetC()) != 0x1a) {
674 return INIT_FAIL_REMOVE;
675 }
676 if ((c = ifs_bitmap->GetC()) == 0x0d) {
677 if ((c = ifs_bitmap->GetC()) != 0x0a) {
678 return INIT_FAIL_REMOVE;
679 }
680 if ((c = ifs_bitmap->GetC()) != 0x1a) {
681 return INIT_FAIL_REMOVE;
682 }
683 if ((c = ifs_bitmap->GetC()) != 0x00) {
684 return INIT_FAIL_REMOVE;
685 }
686 }
687
688 else if (c != 0x00) {
689 return INIT_FAIL_REMOVE;
690 }
691
692 // Read the Color table bit size
693 nColorSize = ifs_bitmap->GetC();
694 if (nColorSize == wxEOF || nColorSize <= 0 || nColorSize > 7) {
695 wxString msg(" Invalid nColorSize data, corrupt on chart ");
696 msg.Append(m_FullPath);
697 wxLogMessage(msg);
698 return INIT_FAIL_REMOVE;
699 }
700
701 // Perform common post-init actions in ChartBaseBSB
702 InitReturn pi_ret = PostInit();
703 if (pi_ret != INIT_OK)
704 return pi_ret;
705 else
706 return INIT_OK;
707}
708
709// ============================================================================
710// ChartKAP implementation
711// ============================================================================
712
713ChartKAP::ChartKAP() { m_ChartType = CHART_TYPE_KAP; }
714
715ChartKAP::~ChartKAP() {}
716
717InitReturn ChartKAP::Init(const wxString &name, ChartInitFlag init_flags) {
718#define BUF_LEN_MAX 4096
719
720 ifs_hdr = new ChartDataNonSeekableInputStream(
721 name); // open the Header file as a read-only stream
722
723 if (!ifs_hdr->IsOk()) return INIT_FAIL_REMOVE;
724
725 int nPlypoint = 0;
726 Plypoint *pPlyTable = (Plypoint *)malloc(sizeof(Plypoint));
727
728 PreInit(name, init_flags, GLOBAL_COLOR_SCHEME_DAY);
729
730 char buffer[BUF_LEN_MAX];
731
732 m_FullPath = name;
733 m_Description = m_FullPath;
734
735 // Clear georeferencing coefficients
736 for (int icl = 0; icl < 12; icl++) {
737 wpx[icl] = 0;
738 wpy[icl] = 0;
739 pwx[icl] = 0;
740 pwy[icl] = 0;
741 }
742
743 // Validate the BSB header
744 // by reading some characters into a buffer and looking for BSB\ keyword
745
746 unsigned int TestBlockSize = 1999;
747 ifs_hdr->Read(buffer, TestBlockSize);
748
749 if (ifs_hdr->LastRead() != TestBlockSize) {
750 wxString msg;
751 msg.Printf(" Could not read first %d bytes of header for chart file: ",
752 TestBlockSize);
753 msg.Append(name);
754 wxLogMessage(msg);
755 free(pPlyTable);
756 return INIT_FAIL_REMOVE;
757 }
758
759 unsigned int i;
760 for (i = 0; i < TestBlockSize - 4; i++) {
761 // Test for "BSB/"
762 if (buffer[i + 0] == 'B' && buffer[i + 1] == 'S' && buffer[i + 2] == 'B' &&
763 buffer[i + 3] == '/')
764 break;
765
766 // Test for "NOS/"
767 if (buffer[i + 0] == 'N' && buffer[i + 1] == 'O' && buffer[i + 2] == 'S' &&
768 buffer[i + 3] == '/')
769 break;
770 }
771 if (i == TestBlockSize - 4) {
772 wxString msg(" Chart file has no BSB header, cannot Init.");
773 msg.Append(name);
774 wxLogMessage(msg);
775 free(pPlyTable);
776 return INIT_FAIL_REMOVE;
777 }
778
779 // Read and Parse Chart Header, line by line
780 ifs_hdr->SeekI(0, wxFromStart); // rewind
781
782 Size_X = Size_Y = 0;
783
784 int done_header_parse = 0;
785 wxCSConv iso_conv("ISO-8859-1"); // we will need a converter
786
787 while (done_header_parse == 0) {
788 if (ReadBSBHdrLine(ifs_hdr, buffer, BUF_LEN_MAX) == 0) {
789 unsigned char c;
790 c = ifs_hdr->GetC();
791 ifs_hdr->Ungetch(c);
792
793 if (0x1a == c)
794 done_header_parse = 1;
795 else {
796 free(pPlyTable);
797 return INIT_FAIL_REMOVE;
798 }
799
800 continue;
801 }
802
803 wxString str_buf(buffer, wxConvUTF8);
804 if (!str_buf.Len()) // failed conversion
805 str_buf = wxString(buffer, iso_conv);
806
807 if (str_buf.Find("SHOM") != wxNOT_FOUND) m_b_SHOM = true;
808
809 if (!strncmp(buffer, "BSB", 3)) {
810 wxString clip_str_buf(
811 &buffer[0],
812 iso_conv); // for single byte French encodings of NAme field
813 wxStringTokenizer tkz(clip_str_buf, "/,=");
814 while (tkz.HasMoreTokens()) {
815 wxString token = tkz.GetNextToken();
816 if (token.IsSameAs("RA", TRUE)) // extract RA=x,y
817 {
818 int i;
819 i = tkz.GetPosition();
820 Size_X = atoi(&buffer[i]);
821 wxString token = tkz.GetNextToken();
822 i = tkz.GetPosition();
823 Size_Y = atoi(&buffer[i]);
824 } else if (token.IsSameAs("NA", TRUE)) // extract NA=str
825 {
826 int i = tkz.GetPosition();
827 char nbuf[81];
828 int j = 0;
829 while ((buffer[i] != ',') && (i < 80)) nbuf[j++] = buffer[i++];
830 nbuf[j] = 0;
831 wxString n_str(nbuf, iso_conv);
832 m_Name = n_str;
833 } else if (token.IsSameAs("NU", TRUE)) // extract NU=str
834 {
835 int i = tkz.GetPosition();
836 char nbuf[81];
837 int j = 0;
838 while ((buffer[i] != ',') && (i < 80)) nbuf[j++] = buffer[i++];
839 nbuf[j] = 0;
840 wxString n_str(nbuf, iso_conv);
841 m_ID = n_str;
842 } else if (token.IsSameAs("DU", TRUE)) // extract DU=n
843 {
844 token = tkz.GetNextToken();
845 long temp_du;
846 if (token.ToLong(&temp_du)) m_Chart_DU = temp_du;
847 }
848 }
849 }
850
851 else if (!strncmp(buffer, "KNP", 3)) {
852 wxString conv_buf(buffer, iso_conv);
853 wxStringTokenizer tkz(conv_buf, "/,=");
854 while (tkz.HasMoreTokens()) {
855 wxString token = tkz.GetNextToken();
856 if (token.IsSameAs("SC", TRUE)) // extract Scale
857 {
858 int i;
859 i = tkz.GetPosition();
860 m_Chart_Scale = atoi(&buffer[i]);
861 if (0 == m_Chart_Scale) m_Chart_Scale = 100000000;
862 } else if (token.IsSameAs("SK", TRUE)) // extract Skew
863 {
864 int i;
865 i = tkz.GetPosition();
866 float fcs;
867 sscanf(&buffer[i], "%f,", &fcs);
868 m_Chart_Skew = fcs;
869 } else if (token.IsSameAs("UN", TRUE)) // extract Depth Units
870 {
871 int i;
872 i = tkz.GetPosition();
873 wxString str(&buffer[i], iso_conv);
874 m_DepthUnits = str.BeforeFirst(',');
875 } else if (token.IsSameAs("GD", TRUE)) // extract Datum
876 {
877 int i;
878 i = tkz.GetPosition();
879 wxString str(&buffer[i], iso_conv);
880 m_datum_str = str.BeforeFirst(',').Trim();
881 } else if (token.IsSameAs("SD", TRUE)) // extract Soundings Datum
882 {
883 int i;
884 i = tkz.GetPosition();
885 wxString str(&buffer[i], iso_conv);
886 m_SoundingsDatum = str.BeforeFirst(',').Trim();
887 } else if (token.IsSameAs("PP",
888 TRUE)) // extract Projection Parameter
889 {
890 int i;
891 i = tkz.GetPosition();
892 double fcs;
893 wxString str(&buffer[i], iso_conv);
894 wxString str1 = str.BeforeFirst(',').Trim();
895 if (str1.ToDouble(&fcs)) m_proj_parameter = fcs;
896 } else if (token.IsSameAs("PR", TRUE)) // extract Projection Type
897 {
898 int i;
899 i = tkz.GetPosition();
900 wxString str(&buffer[i], iso_conv);
901 wxString stru = str.MakeUpper();
902 bool bp_set = false;
903 ;
904
905 if (stru.Matches("*MERCATOR*")) {
906 m_projection = PROJECTION_MERCATOR;
907 bp_set = true;
908 }
909
910 if (stru.Matches("*TRANSVERSE*")) {
911 m_projection = PROJECTION_TRANSVERSE_MERCATOR;
912 bp_set = true;
913 }
914
915 if (stru.Matches("*CONIC*")) {
916 m_projection = PROJECTION_POLYCONIC;
917 bp_set = true;
918 }
919
920 if (stru.Matches("*TM*")) {
921 m_projection = PROJECTION_TRANSVERSE_MERCATOR;
922 bp_set = true;
923 }
924
925 if (stru.Matches("*GAUSS CONFORMAL*")) {
926 m_projection = PROJECTION_TRANSVERSE_MERCATOR;
927 bp_set = true;
928 }
929
930 if (!bp_set) {
931 m_projection = PROJECTION_UNKNOWN;
932 wxString msg(" Chart projection is ");
933 msg += tkz.GetNextToken();
934 msg += " which is unsupported. Disabling chart ";
935 msg += m_FullPath;
936 wxLogMessage(msg);
937 free(pPlyTable);
938 return INIT_FAIL_REMOVE;
939 }
940 } else if (token.IsSameAs(
941 "DX",
942 TRUE)) // extract Pixel scale parameter, if present
943 {
944 int i;
945 i = tkz.GetPosition();
946 float x;
947 sscanf(&buffer[i], "%f,", &x);
948 m_dx = x;
949 } else if (token.IsSameAs(
950 "DY",
951 TRUE)) // extract Pixel scale parameter, if present
952 {
953 int i;
954 i = tkz.GetPosition();
955 float x;
956 sscanf(&buffer[i], "%f,", &x);
957 m_dy = x;
958 }
959 }
960 }
961
962 else if (!strncmp(buffer, "RGB", 3))
963 CreatePaletteEntry(buffer, COLOR_RGB_DEFAULT);
964
965 else if (!strncmp(buffer, "DAY", 3))
966 CreatePaletteEntry(buffer, DAY);
967
968 else if (!strncmp(buffer, "DSK", 3))
969 CreatePaletteEntry(buffer, DUSK);
970
971 else if (!strncmp(buffer, "NGT", 3))
972 CreatePaletteEntry(buffer, NIGHT);
973
974 else if (!strncmp(buffer, "NGR", 3))
975 CreatePaletteEntry(buffer, NIGHTRED);
976
977 else if (!strncmp(buffer, "GRY", 3))
978 CreatePaletteEntry(buffer, GRAY);
979
980 else if (!strncmp(buffer, "PRC", 3))
981 CreatePaletteEntry(buffer, PRC);
982
983 else if (!strncmp(buffer, "PRG", 3))
984 CreatePaletteEntry(buffer, PRG);
985
986 else if (!strncmp(buffer, "REF", 3)) {
987 int i, xr, yr;
988 float ltr, lnr;
989 sscanf(&buffer[4], "%d,%d,%d,%f,%f", &i, &xr, &yr, &ltr, &lnr);
990 pRefTable =
991 (Refpoint *)realloc(pRefTable, sizeof(Refpoint) * (nRefpoint + 1));
992 pRefTable[nRefpoint].xr = xr;
993 pRefTable[nRefpoint].yr = yr;
994 pRefTable[nRefpoint].latr = ltr;
995 pRefTable[nRefpoint].lonr = lnr;
996 pRefTable[nRefpoint].bXValid = 1;
997 pRefTable[nRefpoint].bYValid = 1;
998
999 nRefpoint++;
1000
1001 }
1002
1003 else if (!strncmp(buffer, "WPX", 3)) {
1004 int idx = 0;
1005 double d;
1006 wxStringTokenizer tkz(str_buf.Mid(4), ",");
1007 wxString token = tkz.GetNextToken();
1008
1009 if (token.ToLong((long int *)&wpx_type)) {
1010 while (tkz.HasMoreTokens() && (idx < 12)) {
1011 token = tkz.GetNextToken();
1012 if (token.ToDouble(&d)) {
1013 wpx[idx] = d;
1014 idx++;
1015 }
1016 }
1017 }
1018 n_wpx = idx;
1019 }
1020
1021 else if (!strncmp(buffer, "WPY", 3)) {
1022 int idx = 0;
1023 double d;
1024 wxStringTokenizer tkz(str_buf.Mid(4), ",");
1025 wxString token = tkz.GetNextToken();
1026
1027 if (token.ToLong((long int *)&wpy_type)) {
1028 while (tkz.HasMoreTokens() && (idx < 12)) {
1029 token = tkz.GetNextToken();
1030 if (token.ToDouble(&d)) {
1031 wpy[idx] = d;
1032 idx++;
1033 }
1034 }
1035 }
1036 n_wpy = idx;
1037 }
1038
1039 else if (!strncmp(buffer, "PWX", 3)) {
1040 int idx = 0;
1041 double d;
1042 wxStringTokenizer tkz(str_buf.Mid(4), ",");
1043 wxString token = tkz.GetNextToken();
1044
1045 if (token.ToLong((long int *)&pwx_type)) {
1046 while (tkz.HasMoreTokens() && (idx < 12)) {
1047 token = tkz.GetNextToken();
1048 if (token.ToDouble(&d)) {
1049 pwx[idx] = d;
1050 idx++;
1051 }
1052 }
1053 }
1054 n_pwx = idx;
1055 }
1056
1057 else if (!strncmp(buffer, "PWY", 3)) {
1058 int idx = 0;
1059 double d;
1060 wxStringTokenizer tkz(str_buf.Mid(4), ",");
1061 wxString token = tkz.GetNextToken();
1062
1063 if (token.ToLong((long int *)&pwy_type)) {
1064 while (tkz.HasMoreTokens() && (idx < 12)) {
1065 token = tkz.GetNextToken();
1066 if (token.ToDouble(&d)) {
1067 pwy[idx] = d;
1068 idx++;
1069 }
1070 }
1071 }
1072 n_pwy = idx;
1073 }
1074
1075 else if (!strncmp(buffer, "CPH", 3)) {
1076 float float_cph;
1077 sscanf(&buffer[4], "%f", &float_cph);
1078 m_cph = float_cph;
1079 }
1080
1081 else if (!strncmp(buffer, "VER", 3)) {
1082 wxStringTokenizer tkz(str_buf, "/,=");
1083 wxString token = tkz.GetNextToken();
1084
1085 m_bsb_ver = tkz.GetNextToken();
1086 }
1087
1088 else if (!strncmp(buffer, "DTM", 3)) {
1089 double val;
1090 wxStringTokenizer tkz(str_buf, "/,=");
1091 wxString token = tkz.GetNextToken();
1092
1093 token = tkz.GetNextToken();
1094 if (token.ToDouble(&val)) m_dtm_lat = val;
1095
1096 token = tkz.GetNextToken();
1097 if (token.ToDouble(&val)) m_dtm_lon = val;
1098
1099 // float fdtmlat, fdtmlon;
1100 // sscanf(&buffer[4], "%f,%f", &fdtmlat, &fdtmlon);
1101 // m_dtm_lat = fdtmlat;
1102 // m_dtm_lon = fdtmlon;
1103 }
1104
1105 else if (!strncmp(buffer, "PLY", 3)) {
1106 int i;
1107 float ltp, lnp;
1108 if (sscanf(&buffer[4], "%d,%f,%f", &i, &ltp, &lnp) != 3) {
1109 free(pPlyTable);
1110 return INIT_FAIL_REMOVE;
1111 }
1112 Plypoint *tmp = pPlyTable;
1113 pPlyTable =
1114 (Plypoint *)realloc(pPlyTable, sizeof(Plypoint) * (nPlypoint + 1));
1115 if (NULL == pPlyTable) {
1116 free(tmp);
1117 tmp = NULL;
1118 } else {
1119 pPlyTable[nPlypoint].ltp = ltp;
1120 pPlyTable[nPlypoint].lnp = lnp;
1121 nPlypoint++;
1122 }
1123 if (NULL == pPlyTable || nPlypoint > 1000000) {
1124 // arbitrary 8MB for pPlyTable
1125 nPlypoint = 0;
1126 break;
1127 }
1128 }
1129
1130 else if (!strncmp(buffer, "CED", 3)) {
1131 wxStringTokenizer tkz(str_buf, "/,=");
1132 while (tkz.HasMoreTokens()) {
1133 wxString token = tkz.GetNextToken();
1134 if (token.IsSameAs("ED", TRUE)) // extract Edition Date
1135 {
1136 int i;
1137 i = tkz.GetPosition();
1138
1139 char date_string[40];
1140 char date_buf[16];
1141 date_string[0] = 0;
1142 date_buf[0] = 0;
1143 sscanf(&buffer[i], "%s\r\n", date_string);
1144 wxString date_wxstr(date_string, wxConvUTF8);
1145
1146 wxDateTime dt;
1147 if (dt.ParseDate(date_wxstr)) // successful parse?
1148 {
1149 int iyear =
1150 dt.GetYear(); // GetYear() fails on W98, DMC compiler, wx2.8.3
1151 // BSB charts typically list publish date as xx/yy/zz
1152 // This our own little version of the Y2K problem.
1153 // Just apply some sensible logic
1154
1155 if (iyear < 50) {
1156 iyear += 2000;
1157 dt.SetYear(iyear);
1158 } else if ((iyear >= 50) && (iyear < 100)) {
1159 iyear += 1900;
1160 dt.SetYear(iyear);
1161 }
1162 assert(iyear <= 9999);
1163 sprintf(date_buf, "%d", iyear);
1164
1165 // Initialize the wxDateTime menber for Edition Date
1166 m_EdDate = dt;
1167 } else {
1168 sscanf(date_string, "%s", date_buf);
1169 m_EdDate.Set(1, wxDateTime::Jan,
1170 2000); // Todo this could be smarter
1171 }
1172
1173 m_PubYear = wxString(date_buf, wxConvUTF8);
1174 } else if (token.IsSameAs("SE", TRUE)) // extract Source Edition
1175 {
1176 int i;
1177 i = tkz.GetPosition();
1178 wxString str(&buffer[i], iso_conv);
1179 m_SE = str.BeforeFirst(',');
1180 }
1181 }
1182 }
1183 } // while
1184
1185 // Some charts improperly encode the DTM parameters.
1186 // Identify them as necessary, for further processing
1187 if (m_b_SHOM && (m_bsb_ver == "1.1")) m_b_apply_dtm = false;
1188
1189 // If imbedded coefficients are found,
1190 // then use the polynomial georeferencing algorithms
1191 if (n_pwx && n_pwy && n_pwx && n_pwy) bHaveEmbeddedGeoref = true;
1192
1193 // Set up the projection point according to the projection parameter
1194 if (m_projection == PROJECTION_MERCATOR)
1195 m_proj_lat = m_proj_parameter;
1196 else if (m_projection == PROJECTION_TRANSVERSE_MERCATOR)
1197 m_proj_lon = m_proj_parameter;
1198 else if (m_projection == PROJECTION_POLYCONIC)
1199 m_proj_lon = m_proj_parameter;
1200
1201 // We have seen improperly coded charts, with non-sense value of PP
1202 // parameter FS#1251 Check and override if necessary
1203 if (m_proj_lat > 82.0 || m_proj_lat < -82.0) m_proj_lat = 0.0;
1204
1205 // Validate some of the header data
1206 if (Size_X <= 0 || Size_Y <= 0) {
1207 free(pPlyTable);
1208 return INIT_FAIL_REMOVE;
1209 }
1210
1211 if (nPlypoint < 3) {
1212 wxString msg(" Chart File contains less than 3 or too many PLY points: ");
1213 msg.Append(m_FullPath);
1214 wxLogMessage(msg);
1215 free(pPlyTable);
1216 return INIT_FAIL_REMOVE;
1217 }
1218
1219 if (m_datum_str.IsEmpty()) {
1220 wxString msg(" Chart datum not specified on chart ");
1221 msg.Append(m_FullPath);
1222 wxLogMessage(msg);
1223 wxLogMessage(" Default datum (WGS84) substituted.");
1224
1225 // return INIT_FAIL_REMOVE;
1226 } else {
1227 char d_str[100];
1228 strncpy(d_str, m_datum_str.mb_str(), 99);
1229 d_str[99] = 0;
1230
1231 int datum_index = GetDatumIndex(d_str);
1232
1233 if (datum_index < 0) {
1234 wxString msg(" Chart datum {");
1235 msg += m_datum_str;
1236 msg += "} invalid on chart ";
1237 msg.Append(m_FullPath);
1238 wxLogMessage(msg);
1239 wxLogMessage(" Default datum (WGS84) substituted.");
1240
1241 // return INIT_FAIL_REMOVE;
1242 }
1243 }
1244
1245 /* Augment ply points
1246 This is needed for example on polyconic charts or skewed charts because
1247 straight lines in the chart coordinates can not use simple
1248 interpolation in lat/lon or mercator coordinate space to draw the
1249 borders or be used for quilting operation.
1250 TODO: should this be added as a subroutine for GEO chartso? */
1251 if ((m_projection != PROJECTION_MERCATOR &&
1252 m_projection != PROJECTION_TRANSVERSE_MERCATOR) ||
1253 m_Chart_Skew > 2) {
1254 // Analyze Refpoints early because we need georef coefficient here.
1255 AnalyzeRefpoints(false); // no post test needed
1256
1257 // We need to compute a tentative min/max lat/lon to perform georefs
1258 // These lat/lon extents will be more accurately updated later.
1259 m_LonMax = -360.0;
1260 m_LonMin = 360.0;
1261 m_LatMax = -90.0;
1262 m_LatMin = 90.0;
1263
1264 for (int i = 0; i < nPlypoint; i++) {
1265 m_LatMax = wxMax(m_LatMax, pPlyTable[i].ltp);
1266 m_LatMin = wxMin(m_LatMin, pPlyTable[i].ltp);
1267 m_LonMax = wxMax(m_LonMax, pPlyTable[i].lnp);
1268 m_LonMin = wxMin(m_LonMin, pPlyTable[i].lnp);
1269 }
1270
1271 int count = nPlypoint;
1272 nPlypoint = 0;
1273 Plypoint *pOldPlyTable = pPlyTable;
1274 pPlyTable = NULL;
1275 double lastplylat = 0.0, lastplylon = 0.0, x1 = 0.0, y1 = 0.0, x2, y2;
1276 double plylat, plylon;
1277 for (int i = 0; i < count + 1; i++) {
1278 plylat = pOldPlyTable[i % count].ltp;
1279 plylon = pOldPlyTable[i % count].lnp;
1280 latlong_to_chartpix(plylat, plylon, x2, y2);
1281 if (i > 0) {
1282 if (lastplylon - plylon > 180.)
1283 lastplylon -= 360.;
1284 else if (lastplylon - plylon < -180.)
1285 lastplylon += 360.;
1286
1287 // use 2 degree steps
1288 double steps =
1289 ceil((fabs(lastplylat - plylat) + fabs(lastplylon - plylon)) / 2);
1290 for (double c = 0; c < steps; c++) {
1291 double d = c / steps, lat, lon;
1292 wxPoint2DDouble s;
1293 double x = (1 - d) * x1 + d * x2, y = (1 - d) * y1 + d * y2;
1294 chartpix_to_latlong(x, y, &lat, &lon);
1295 pPlyTable = (Plypoint *)realloc(pPlyTable,
1296 sizeof(Plypoint) * (nPlypoint + 1));
1297 pPlyTable[nPlypoint].ltp = lat;
1298 pPlyTable[nPlypoint].lnp = lon;
1299 nPlypoint++;
1300 }
1301 }
1302 x1 = x2, y1 = y2;
1303 lastplylat = plylat, lastplylon = plylon;
1304 }
1305 free(pOldPlyTable);
1306 }
1307
1308 // Convert captured plypoint information into chart COVR structures
1309
1310 // A special-case test for poorly formatted charts
1311 // We look for cases where the declared PlyPoints are far outside of the
1312 // chart raster bitmap If found, we change the COVR region to the valid
1313 // bitmap region, instead of the default PlyPoints region
1314 // Set a tentative lat/lon range.
1315 m_LonMax = -360.;
1316 m_LonMin = 360.;
1317 for (int i = 0; i < nPlypoint; i++) {
1318 m_LonMin = wxMin(m_LonMin, pPlyTable[i].lnp);
1319 m_LonMax = wxMax(m_LonMax, pPlyTable[i].lnp);
1320 }
1321 // This test does not really work for charts that cross IDL
1322 bool b_test = true;
1323 bool b_adjusted = false;
1324 if (m_LonMax * m_LonMin < 0) {
1325 if ((m_LonMax - m_LonMin) > 180.) b_test = false;
1326 }
1327
1328 if (b_test) {
1329 if (!bHaveEmbeddedGeoref) {
1330 // Analyze Refpoints early because we might need georef coefficient
1331 // here.
1332 AnalyzeRefpoints(false); // no post test needed
1333 }
1334
1335 bool bAdjustPly = false;
1336 wxRect bitRect(0, 0, Size_X, Size_Y);
1337 bitRect.Inflate(5); // allow for a little roundoff error
1338 for (int i = 0; i < nPlypoint; i++) {
1339 double pix_x, pix_y;
1340 latlong_to_chartpix(pPlyTable[i].ltp, pPlyTable[i].lnp, pix_x, pix_y);
1341 if (!bitRect.Contains(pix_x, pix_y)) {
1342 bAdjustPly = true;
1343 if (m_b_cdebug)
1344 printf("Adjusting COVR region on: %s\n", name.ToUTF8().data());
1345 break;
1346 }
1347 }
1348
1349 if (bAdjustPly) {
1350 float *points = new float[2 * nPlypoint];
1351 for (int i = 0; i < nPlypoint; i++)
1352 points[2 * i + 0] = pPlyTable[i].ltp,
1353 points[2 * i + 1] = pPlyTable[i].lnp;
1354 LLRegion covrRegion(nPlypoint, points);
1355 delete[] points;
1356 covrRegion.Intersect(GetValidRegion());
1357
1358 if (covrRegion.contours.size()) { // Check for no intersection caused by
1359 // bogus georef....
1360 m_nCOVREntries = covrRegion.contours.size();
1361 m_pCOVRTablePoints = (int *)malloc(m_nCOVREntries * sizeof(int));
1362 m_pCOVRTable = (float **)malloc(m_nCOVREntries * sizeof(float *));
1363 std::list<poly_contour>::iterator it = covrRegion.contours.begin();
1364 for (int i = 0; i < m_nCOVREntries; i++) {
1365 m_pCOVRTablePoints[i] = it->size();
1366 m_pCOVRTable[i] =
1367 (float *)malloc(m_pCOVRTablePoints[i] * 2 * sizeof(float));
1368 std::list<contour_pt>::iterator jt = it->begin();
1369 for (int j = 0; j < m_pCOVRTablePoints[i]; j++) {
1370 m_pCOVRTable[i][2 * j + 0] = jt->y;
1371 m_pCOVRTable[i][2 * j + 1] = jt->x;
1372 jt++;
1373 }
1374 it++;
1375 }
1376 b_adjusted = true;
1377 }
1378 }
1379 }
1380
1381 if (!b_adjusted) {
1382 m_nCOVREntries = 1;
1383 m_pCOVRTablePoints = (int *)malloc(sizeof(int));
1384 *m_pCOVRTablePoints = nPlypoint;
1385 m_pCOVRTable = (float **)malloc(sizeof(float *));
1386 *m_pCOVRTable = (float *)malloc(nPlypoint * 2 * sizeof(float));
1387 memcpy(*m_pCOVRTable, pPlyTable, nPlypoint * 2 * sizeof(float));
1388 }
1389
1390 free(pPlyTable);
1391
1392 // Setup the datum transform parameters
1393 char d_str[100];
1394 strncpy(d_str, m_datum_str.mb_str(), 99);
1395 d_str[99] = 0;
1396
1397 int datum_index = GetDatumIndex(d_str);
1398 m_datum_index = datum_index;
1399
1400 if (datum_index < 0)
1401 m_ExtraInfo = "---<<< Warning: Chart Datum may be incorrect. >>>---";
1402
1403 // Establish defaults, may be overridden later
1404 m_lon_datum_adjust = (-m_dtm_lon) / 3600.;
1405 m_lat_datum_adjust = (-m_dtm_lat) / 3600.;
1406
1407 // Adjust the PLY points to WGS84 datum
1408 Plypoint *ppp = (Plypoint *)GetCOVRTableHead(0);
1409 int cnPlypoint = GetCOVRTablenPoints(0);
1410
1411 for (int u = 0; u < cnPlypoint; u++) {
1412 double dlat = 0;
1413 double dlon = 0;
1414
1415 if (m_datum_index == DATUM_INDEX_WGS84 ||
1416 m_datum_index == DATUM_INDEX_UNKNOWN) {
1417 dlon = m_dtm_lon / 3600.;
1418 dlat = m_dtm_lat / 3600.;
1419 }
1420
1421 else {
1422 double to_lat, to_lon;
1423 MolodenskyTransform(ppp->ltp, ppp->lnp, &to_lat, &to_lon, m_datum_index,
1424 DATUM_INDEX_WGS84);
1425 dlon = (to_lon - ppp->lnp);
1426 dlat = (to_lat - ppp->ltp);
1427 if (m_b_apply_dtm) {
1428 dlon += m_dtm_lon / 3600.;
1429 dlat += m_dtm_lat / 3600.;
1430 }
1431 }
1432
1433 ppp->lnp += dlon;
1434 ppp->ltp += dlat;
1435 ppp++;
1436 }
1437
1438 if (!SetMinMax()) return INIT_FAIL_REMOVE; // have to bail here
1439
1440 AnalyzeSkew();
1441
1442 if (init_flags == HEADER_ONLY) return INIT_OK;
1443
1444 // Advance to the data
1445 unsigned char c;
1446 bool bcorrupt = false;
1447
1448 if ((c = ifs_hdr->GetC()) != 0x1a) {
1449 bcorrupt = true;
1450 }
1451 if ((c = ifs_hdr->GetC()) == 0x0d) {
1452 if ((c = ifs_hdr->GetC()) != 0x0a) {
1453 bcorrupt = true;
1454 }
1455 if ((c = ifs_hdr->GetC()) != 0x1a) {
1456 bcorrupt = true;
1457 }
1458 if ((c = ifs_hdr->GetC()) != 0x00) {
1459 bcorrupt = true;
1460 }
1461 }
1462
1463 else if (c != 0x00) {
1464 bcorrupt = true;
1465 }
1466
1467 if (bcorrupt) {
1468 wxString msg(" Chart File RLL data corrupt on chart ");
1469 msg.Append(m_FullPath);
1470 wxLogMessage(msg);
1471
1472 return INIT_FAIL_REMOVE;
1473 }
1474
1475 // Read the Color table bit size
1476 nColorSize = ifs_hdr->GetC();
1477 if (nColorSize == wxEOF || nColorSize <= 0 || nColorSize > 7) {
1478 wxString msg(" Invalid nColorSize data, corrupt on chart ");
1479 msg.Append(m_FullPath);
1480 wxLogMessage(msg);
1481 return INIT_FAIL_REMOVE;
1482 }
1483
1484 nFileOffsetDataStart = ifs_hdr->TellI();
1485 delete ifs_hdr;
1486 ifs_hdr = NULL;
1487
1488 ChartDataInputStream *stream =
1489 new ChartDataInputStream(name); // Open again, as the bitmap
1490 wxString tempfile;
1491#ifdef OCPN_USE_LZMA
1492 tempfile = stream->TempFileName();
1493#endif
1494 m_filesize = wxFileName::GetSize(tempfile.empty() ? name : tempfile);
1495
1496 ifss_bitmap = stream;
1497 ifs_bitmap = new wxBufferedInputStream(*ifss_bitmap);
1498
1499 // Perform common post-init actions in ChartBaseBSB
1500 InitReturn pi_ret = PostInit();
1501 if (pi_ret != INIT_OK) return pi_ret;
1502 return INIT_OK;
1503}
1504
1505// ============================================================================
1506// ChartBaseBSB implementation
1507// ============================================================================
1508
1509ChartBaseBSB::ChartBaseBSB() {
1510 // Init some private data
1511 m_ChartFamily = CHART_FAMILY_RASTER;
1512
1513 pBitmapFilePath = NULL;
1514
1515 pline_table = NULL;
1516 ifs_buf = NULL;
1517
1518 cached_image_ok = 0;
1519
1520 pRefTable = (Refpoint *)malloc(sizeof(Refpoint));
1521 nRefpoint = 0;
1522 cPoints.status = 0;
1523 bHaveEmbeddedGeoref = false;
1524 n_wpx = 0;
1525 n_wpy = 0;
1526 n_pwx = 0;
1527 n_pwy = 0;
1528
1529#ifdef __ANDROID__
1530 bUseLineCache = false;
1531#else
1532 bUseLineCache = true;
1533#endif
1534
1535 m_Chart_Skew = 0.0;
1536
1537 pPixCache = NULL;
1538
1539 pLineCache = NULL;
1540
1541 m_bilinear_limit = 8; // bilinear scaling only up to n
1542
1543 ifs_bitmap = NULL;
1544 ifss_bitmap = NULL;
1545 ifs_hdr = NULL;
1546
1547 for (int i = 0; i < N_BSB_COLORS; i++) pPalettes[i] = NULL;
1548
1549 bGeoErrorSent = false;
1550 m_Chart_DU = 0;
1551 m_cph = 0.;
1552
1553 m_mapped_color_index = COLOR_RGB_DEFAULT;
1554
1555 m_datum_str = "WGS84"; // assume until proven otherwise
1556
1557 m_dtm_lat = 0.;
1558 m_dtm_lon = 0.;
1559
1560 m_dx = 0.;
1561 m_dy = 0.;
1562 m_proj_lat = 0.;
1563 m_proj_lon = 0.;
1564 m_proj_parameter = 0.;
1565 m_b_SHOM = false;
1566 m_b_apply_dtm = true;
1567
1568 m_b_cdebug = 0;
1569
1570#ifdef OCPN_USE_CONFIG
1571 wxFileConfig *pfc = (wxFileConfig *)pConfig;
1572 pfc->SetPath(_T ( "/Settings" ));
1573 pfc->Read(_T ( "DebugBSBImg" ), &m_b_cdebug, 0);
1574#endif
1575}
1576
1577ChartBaseBSB::~ChartBaseBSB() {
1578 if (pBitmapFilePath) delete pBitmapFilePath;
1579
1580 if (pline_table) free(pline_table);
1581
1582 if (ifs_buf) free(ifs_buf);
1583
1584 free(pRefTable);
1585 // free(pPlyTable);
1586
1587 delete ifs_bitmap;
1588 delete ifs_hdr;
1589 delete ifss_bitmap;
1590
1591 if (cPoints.status) {
1592 free(cPoints.tx);
1593 free(cPoints.ty);
1594 free(cPoints.lon);
1595 free(cPoints.lat);
1596
1597 free(cPoints.pwx);
1598 free(cPoints.wpx);
1599 free(cPoints.pwy);
1600 free(cPoints.wpy);
1601 }
1602
1603 // Free the line cache
1604 FreeLineCacheRows();
1605 free(pLineCache);
1606
1607 delete pPixCache;
1608
1609 for (int i = 0; i < N_BSB_COLORS; i++) delete pPalettes[i];
1610}
1611
1612void ChartBaseBSB::FreeLineCacheRows(int start, int end) {
1613 if (pLineCache) {
1614 if (end < 0)
1615 end = Size_Y;
1616 else
1617 end = wxMin(end, Size_Y);
1618 for (int ylc = start; ylc < end; ylc++) {
1619 CachedLine *pt = &pLineCache[ylc];
1620 if (pt->bValid) {
1621 free(pt->pTileOffset);
1622 free(pt->pPix);
1623 pt->bValid = false;
1624 }
1625 }
1626 }
1627}
1628
1629bool ChartBaseBSB::HaveLineCacheRow(int row) {
1630 if (pLineCache) {
1631 CachedLine *pt = &pLineCache[row];
1632 return pt->bValid;
1633 }
1634 return false;
1635}
1636
1637// Report recommended minimum and maximum scale values for which use of this
1638// chart is valid
1639
1640double ChartBaseBSB::GetNormalScaleMin(double canvas_scale_factor,
1641 bool b_allow_overzoom) {
1642 // if(b_allow_overzoom)
1643 return (canvas_scale_factor / m_ppm_avg) /
1644 32; // allow wide range overzoom overscale
1645 // else
1646 // return (canvas_scale_factor / m_ppm_avg) / 2; // don't
1647 // suggest too much overscale
1648}
1649
1650double ChartBaseBSB::GetNormalScaleMax(double canvas_scale_factor,
1651 int canvas_width) {
1652 return (canvas_scale_factor / m_ppm_avg) *
1653 4.0; // excessive underscale is slow, and unreadable
1654}
1655
1656double ChartBaseBSB::GetNearestPreferredScalePPM(double target_scale_ppm) {
1658 target_scale_ppm, .01, 64.); // changed from 32 to 64 to allow super
1659 // small scale BSB charts as quilt base
1660}
1661
1663 double scale_factor_min,
1664 double scale_factor_max) {
1665 double chart_1x_scale = GetPPM();
1666
1667 double binary_scale_factor = 1.;
1668
1669 // Overzoom....
1670 if (chart_1x_scale > target_scale) {
1671 double binary_scale_factor_max = 1 / scale_factor_min;
1672
1673 while (binary_scale_factor < binary_scale_factor_max) {
1674 if (fabs((chart_1x_scale / binary_scale_factor) - target_scale) <
1675 (target_scale * 0.05))
1676 break;
1677 if ((chart_1x_scale / binary_scale_factor) < target_scale)
1678 break;
1679 else
1680 binary_scale_factor *= 2.;
1681 }
1682 }
1683
1684 // Underzoom.....
1685 else {
1686 int ibsf = 1;
1687 int isf_max = (int)scale_factor_max;
1688 while (ibsf < isf_max) {
1689 if (fabs((chart_1x_scale * ibsf) - target_scale) < (target_scale * 0.05))
1690 break;
1691
1692 else if ((chart_1x_scale * ibsf) > target_scale) {
1693 if (ibsf > 1) ibsf /= 2;
1694 break;
1695 } else
1696 ibsf *= 2;
1697 }
1698
1699 binary_scale_factor = 1. / ibsf;
1700 }
1701
1702 return chart_1x_scale / binary_scale_factor;
1703}
1704
1705InitReturn ChartBaseBSB::Init(const wxString &name, ChartInitFlag init_flags) {
1706 m_global_color_scheme = GLOBAL_COLOR_SCHEME_RGB;
1707 return INIT_OK;
1708}
1709
1710InitReturn ChartBaseBSB::PreInit(const wxString &name, ChartInitFlag init_flags,
1711 ColorScheme cs) {
1712 m_global_color_scheme = cs;
1713 return INIT_OK;
1714}
1715
1716void ChartBaseBSB::CreatePaletteEntry(char *buffer, int palette_index) {
1717 if (palette_index < N_BSB_COLORS) {
1718 if (!pPalettes[palette_index]) pPalettes[palette_index] = new opncpnPalette;
1719 opncpnPalette *pp = pPalettes[palette_index];
1720
1721 pp->FwdPalette =
1722 (int *)realloc(pp->FwdPalette, (pp->nFwd + 1) * sizeof(int));
1723 pp->RevPalette =
1724 (int *)realloc(pp->RevPalette, (pp->nRev + 1) * sizeof(int));
1725 pp->nFwd++;
1726 pp->nRev++;
1727
1728 int i;
1729 int n, r, g, b;
1730 sscanf(&buffer[4], "%d,%d,%d,%d", &n, &r, &g, &b);
1731
1732 i = n;
1733
1734 int fcolor, rcolor;
1735 fcolor = (b << 16) + (g << 8) + r;
1736 rcolor = (r << 16) + (g << 8) + b;
1737
1738 pp->RevPalette[i] = rcolor;
1739 pp->FwdPalette[i] = fcolor;
1740 }
1741}
1742
1743InitReturn ChartBaseBSB::PostInit(void) {
1744 // catch undefined shift if not already done in derived classes
1745 if (nColorSize == wxEOF || nColorSize <= 0 || nColorSize > 7) {
1746 wxString msg(" Invalid nColorSize data, corrupt in PostInit() on chart ");
1747 msg.Append(m_FullPath);
1748 wxLogMessage(msg);
1749 return INIT_FAIL_REMOVE;
1750 }
1751
1752 if (Size_X <= 0 || Size_X > INT_MAX / 4 || Size_Y <= 0 ||
1753 Size_Y - 1 > INT_MAX / 4) {
1754 wxString msg(
1755 " Invalid Size_X/Size_Y data, corrupt in PostInit() on chart ");
1756 msg.Append(m_FullPath);
1757 wxLogMessage(msg);
1758 return INIT_FAIL_REMOVE;
1759 }
1760
1761 // Validate the palette array, substituting DEFAULT for missing entries
1762 int nfwd_def = 1;
1763 int nrev_def = 1;
1764 if (pPalettes[COLOR_RGB_DEFAULT]) {
1765 nrev_def = pPalettes[COLOR_RGB_DEFAULT]->nRev;
1766 nfwd_def = pPalettes[COLOR_RGB_DEFAULT]->nFwd;
1767 }
1768
1769 for (int i = 0; i < N_BSB_COLORS; i++) {
1770 if (pPalettes[i] == NULL) {
1771 opncpnPalette *pNullSubPal = new opncpnPalette;
1772
1773 pNullSubPal->nFwd = nfwd_def; // copy the palette count
1774 pNullSubPal->nRev = nrev_def; // copy the palette count
1775 // Deep copy the palette rgb tables
1776 free(pNullSubPal->FwdPalette);
1777 pNullSubPal->FwdPalette = (int *)malloc(pNullSubPal->nFwd * sizeof(int));
1778 if (pPalettes[COLOR_RGB_DEFAULT])
1779 memcpy(pNullSubPal->FwdPalette,
1780 pPalettes[COLOR_RGB_DEFAULT]->FwdPalette,
1781 pNullSubPal->nFwd * sizeof(int));
1782
1783 free(pNullSubPal->RevPalette);
1784 pNullSubPal->RevPalette = (int *)malloc(pNullSubPal->nRev * sizeof(int));
1785 if (pPalettes[COLOR_RGB_DEFAULT])
1786 memcpy(pNullSubPal->RevPalette,
1787 pPalettes[COLOR_RGB_DEFAULT]->RevPalette,
1788 pNullSubPal->nRev * sizeof(int));
1789
1790 pPalettes[i] = pNullSubPal;
1791 }
1792 }
1793
1794 // Establish the palette type and default palette
1795 palette_direction = GetPaletteDir();
1796
1797 SetColorScheme(m_global_color_scheme, false);
1798
1799 // Allocate memory for ifs file buffering
1800 ifs_bufsize = Size_X * 4;
1801 ifs_buf = (unsigned char *)malloc(ifs_bufsize);
1802 if (!ifs_buf) return INIT_FAIL_REMOVE;
1803
1804 ifs_bufend = ifs_buf + ifs_bufsize;
1805 ifs_lp = ifs_bufend;
1806 ifs_file_offset = -ifs_bufsize;
1807
1808 // Create and load the line offset index table
1809 pline_table = NULL;
1810 pline_table = (int *)malloc((Size_Y + 1) * sizeof(int)); // Ugly....
1811 if (!pline_table) return INIT_FAIL_REMOVE;
1812
1813 ifs_bitmap->SeekI((Size_Y + 1) * -4,
1814 wxFromEnd); // go to Beginning of offset table
1815 pline_table[Size_Y] = ifs_bitmap->TellI(); // fill in useful last table entry
1816
1817 unsigned char *tmp = (unsigned char *)malloc(Size_Y * sizeof(int));
1818 ifs_bitmap->Read(tmp, Size_Y * sizeof(int));
1819 if (ifs_bitmap->LastRead() != Size_Y * sizeof(int)) {
1820 wxString msg(" Chart File corrupt in PostInit() on chart ");
1821 msg.Append(m_FullPath);
1822 wxLogMessage(msg);
1823 free(tmp);
1824
1825 return INIT_FAIL_REMOVE;
1826 }
1827
1828 int offset;
1829 unsigned char *b = tmp;
1830 for (int ifplt = 0; ifplt < Size_Y; ifplt++) {
1831 offset = 0;
1832 offset += *b++ * 256 * 256 * 256;
1833 offset += *b++ * 256 * 256;
1834 offset += *b++ * 256;
1835 offset += *b++;
1836
1837 pline_table[ifplt] = offset;
1838 }
1839 free(tmp);
1840 // Try to validate the line index
1841
1842 bool bline_index_ok = true;
1843 m_nLineOffset = 0;
1844
1845 wxULongLong bitmap_filesize = m_filesize;
1846 if ((m_ChartType == CHART_TYPE_GEO) && pBitmapFilePath)
1847 bitmap_filesize = wxFileName::GetSize(*pBitmapFilePath);
1848
1849 // look logically at the line offset table
1850 for (int iplt = 0; iplt < Size_Y - 1; iplt++) {
1851 if (pline_table[iplt] > bitmap_filesize) {
1852 wxString msg(" Chart File corrupt in PostInit() on chart ");
1853 msg.Append(m_FullPath);
1854 wxLogMessage(msg);
1855
1856 return INIT_FAIL_REMOVE;
1857 }
1858
1859 int thisline_size = pline_table[iplt + 1] - pline_table[iplt];
1860 if (thisline_size < 0) {
1861 wxString msg(" Chart File corrupt in PostInit() on chart ");
1862 msg.Append(m_FullPath);
1863 wxLogMessage(msg);
1864
1865 return INIT_FAIL_REMOVE;
1866 }
1867 }
1868
1869 // For older charts, say Version 1.x, we will try to read the chart and check
1870 // the lines for coherence These older charts are more likely to have index
1871 // troubles.... We only need to check a few lines. Errors are quickly
1872 // apparent.
1873 double ver;
1874 m_bsb_ver.ToDouble(&ver);
1875 if (ver < 2.0) {
1876 for (int iplt = 0; iplt < 10; iplt++) {
1877 if (wxInvalidOffset ==
1878 ifs_bitmap->SeekI(pline_table[iplt], wxFromStart)) {
1879 wxString msg(" Chart File corrupt in PostInit() on chart ");
1880 msg.Append(m_FullPath);
1881 wxLogMessage(msg);
1882
1883 return INIT_FAIL_REMOVE;
1884 }
1885
1886 int thisline_size = pline_table[iplt + 1] - pline_table[iplt];
1887 ifs_bitmap->Read(ifs_buf, thisline_size);
1888
1889 unsigned char *lp = ifs_buf;
1890
1891 unsigned char byNext;
1892 int nLineMarker = 0;
1893 do {
1894 byNext = *lp++;
1895 nLineMarker = nLineMarker * 128 + (byNext & 0x7f);
1896 } while ((byNext & 0x80) != 0);
1897
1898 // Linemarker Correction factor needed here
1899 // Some charts start with LineMarker = 0, some with LineMarker = 1
1900 // Assume the first LineMarker found is the index base, and use
1901 // as a correction offset
1902
1903 if (iplt == 0) m_nLineOffset = nLineMarker;
1904
1905 if (nLineMarker != iplt + m_nLineOffset) {
1906 bline_index_ok = false;
1907 break;
1908 }
1909 }
1910 }
1911
1912 // Recreate the scan line index if the embedded version seems corrupt
1913 if (!bline_index_ok) {
1914 wxString msg(" Line Index corrupt, recreating Index for chart ");
1915 msg.Append(m_FullPath);
1916 wxLogMessage(msg);
1917 if (!CreateLineIndex()) {
1918 wxString msg(" Error creating Line Index for chart ");
1919 msg.Append(m_FullPath);
1920 wxLogMessage(msg);
1921 return INIT_FAIL_REMOVE;
1922 }
1923 }
1924
1925 // Allocate the Line Cache
1926 if (bUseLineCache) {
1927 pLineCache = (CachedLine *)malloc(Size_Y * sizeof(CachedLine));
1928 CachedLine *pt;
1929
1930 for (int ylc = 0; ylc < Size_Y; ylc++) {
1931 pt = &pLineCache[ylc];
1932 pt->bValid = false;
1933 pt->pPix = NULL; //(unsigned char *)malloc(1);
1934 pt->pTileOffset = NULL;
1935 }
1936 } else
1937 pLineCache = NULL;
1938
1939 // Validate/Set Depth Unit Type
1940 wxString test_str = m_DepthUnits.Upper();
1941 if (test_str.IsSameAs("FEET", FALSE))
1942 m_depth_unit_id = DEPTH_UNIT_FEET;
1943 else if (test_str.IsSameAs("METERS", FALSE))
1944 m_depth_unit_id = DEPTH_UNIT_METERS;
1945 else if (test_str.IsSameAs("METRES",
1946 FALSE)) // Special case for alternate spelling
1947 m_depth_unit_id = DEPTH_UNIT_METERS;
1948 else if (test_str.IsSameAs("METRIC", FALSE))
1949 m_depth_unit_id = DEPTH_UNIT_METERS;
1950 else if (test_str.IsSameAs("FATHOMS", FALSE))
1951 m_depth_unit_id = DEPTH_UNIT_FATHOMS;
1952 else if (test_str.Find("FATHOMS") !=
1953 wxNOT_FOUND) // Special case for "Fathoms and Feet"
1954 m_depth_unit_id = DEPTH_UNIT_FATHOMS;
1955 else if (test_str.Find("METERS") !=
1956 wxNOT_FOUND) // Special case for "Meters and decimeters"
1957 m_depth_unit_id = DEPTH_UNIT_METERS;
1958
1959 // Analyze Refpoints
1960 int analyze_ret_val = AnalyzeRefpoints();
1961 if (0 != analyze_ret_val) return INIT_FAIL_REMOVE;
1962
1963 bReadyToRender = true;
1964 return INIT_OK;
1965}
1966
1967bool ChartBaseBSB::CreateLineIndex() {
1968 // Assumes file stream ifs_bitmap is currently open
1969
1970 // wxBufferedInputStream *pbis = new wxBufferedInputStream(*ifss_bitmap);
1971
1972 // Seek to start of data
1973 ifs_bitmap->SeekI(nFileOffsetDataStart); // go to Beginning of data
1974
1975 for (int iplt = 0; iplt < Size_Y; iplt++) {
1976 int offset = ifs_bitmap->TellI();
1977
1978 int iscan = BSBScanScanline(ifs_bitmap);
1979
1980 // There is no sense reporting an error here, since we are recreating after
1981 // an error
1982 /*
1983 if(iscan > Size_Y)
1984 {
1985
1986 wxString msg("CreateLineIndex() failed on chart ");
1987 msg.Append(m_FullPath);
1988 wxLogMessage(msg);
1989 return false;
1990 }
1991
1992 // Skipped lines?
1993 if(iscan != iplt)
1994 {
1995 while((iplt < iscan) && (iplt < Size_Y))
1996 {
1997 pline_table[iplt] = 0;
1998 iplt++;
1999 }
2000 }
2001 */
2002 pline_table[iplt] = offset;
2003 }
2004
2005 return true;
2006}
2007
2008// Invalidate and Free the line cache contents
2009void ChartBaseBSB::InvalidateLineCache(void) {
2010 if (pLineCache) {
2011 CachedLine *pt;
2012 for (int ylc = 0; ylc < Size_Y; ylc++) {
2013 pt = &pLineCache[ylc];
2014 if (pt) {
2015 free(pt->pPix);
2016 pt->pPix = NULL;
2017 free(pt->pTileOffset);
2018 pt->pTileOffset = NULL;
2019 pt->bValid = false;
2020 }
2021 }
2022 }
2023}
2024
2025bool ChartBaseBSB::GetChartExtent(Extent *pext) {
2026 pext->NLAT = m_LatMax;
2027 pext->SLAT = m_LatMin;
2028 pext->ELON = m_LonMax;
2029 pext->WLON = m_LonMin;
2030
2031 return true;
2032}
2033
2034bool ChartBaseBSB::SetMinMax(void) {
2035 // Calculate the Chart Extents(M_LatMin, M_LonMin, etc.)
2036 // from the COVR data, for fast database search
2037 m_LonMax = -360.0;
2038 m_LonMin = 360.0;
2039 m_LatMax = -90.0;
2040 m_LatMin = 90.0;
2041
2042 Plypoint *ppp = (Plypoint *)GetCOVRTableHead(0);
2043 int cnPlypoint = GetCOVRTablenPoints(0);
2044
2045 for (int u = 0; u < cnPlypoint; u++) {
2046 if (ppp->lnp > m_LonMax) m_LonMax = ppp->lnp;
2047 if (ppp->lnp < m_LonMin) m_LonMin = ppp->lnp;
2048
2049 if (ppp->ltp > m_LatMax) m_LatMax = ppp->ltp;
2050 if (ppp->ltp < m_LatMin) m_LatMin = ppp->ltp;
2051
2052 ppp++;
2053 }
2054
2055 // Check for special cases
2056
2057 // Case 1: Chart spans International Date Line or Greenwich, Longitude
2058 // min/max is non-obvious.
2059 if ((m_LonMax * m_LonMin) < 0) // min/max are opposite signs
2060 {
2061 // Georeferencing is not yet available, so find the reference points
2062 // closest to min/max ply points
2063
2064 if (0 == nRefpoint) return false; // have to bail here
2065
2066 // for m_LonMax
2067 double min_dist_x = 360;
2068 int imaxclose = 0;
2069 for (int ic = 0; ic < nRefpoint; ic++) {
2070 double dist = sqrt(
2071 ((m_LatMax - pRefTable[ic].latr) * (m_LatMax - pRefTable[ic].latr)) +
2072 ((m_LonMax - pRefTable[ic].lonr) * (m_LonMax - pRefTable[ic].lonr)));
2073
2074 if (dist < min_dist_x) {
2075 min_dist_x = dist;
2076 imaxclose = ic;
2077 }
2078 }
2079
2080 // for m_LonMin
2081 double min_dist_n = 360;
2082 int iminclose = 0;
2083 for (int id = 0; id < nRefpoint; id++) {
2084 double dist = sqrt(
2085 ((m_LatMin - pRefTable[id].latr) * (m_LatMin - pRefTable[id].latr)) +
2086 ((m_LonMin - pRefTable[id].lonr) * (m_LonMin - pRefTable[id].lonr)));
2087
2088 if (dist < min_dist_n) {
2089 min_dist_n = dist;
2090 iminclose = id;
2091 }
2092 }
2093
2094 // Is this chart crossing IDL or Greenwich?
2095 // Make the check
2096 if (pRefTable[imaxclose].xr < pRefTable[iminclose].xr) {
2097 // This chart crosses IDL and needs a flip, meaning that all negative
2098 // longitudes need to be normalized and the min/max relcalculated This
2099 // code added to correct non-rectangular charts crossing IDL, such as
2100 // nz14605.kap
2101
2102 m_LonMax = -360.0;
2103 m_LonMin = 360.0;
2104 m_LatMax = -90.0;
2105 m_LatMin = 90.0;
2106
2107 Plypoint *ppp =
2108 (Plypoint *)GetCOVRTableHead(0); // Normalize the plypoints
2109 int cnPlypoint = GetCOVRTablenPoints(0);
2110
2111 for (int u = 0; u < cnPlypoint; u++) {
2112 if (ppp->lnp < 0.) ppp->lnp += 360.;
2113
2114 if (ppp->lnp > m_LonMax) m_LonMax = ppp->lnp;
2115 if (ppp->lnp < m_LonMin) m_LonMin = ppp->lnp;
2116
2117 if (ppp->ltp > m_LatMax) m_LatMax = ppp->ltp;
2118 if (ppp->ltp < m_LatMin) m_LatMin = ppp->ltp;
2119
2120 ppp++;
2121 }
2122 }
2123 }
2124
2125 // Case 2 Lons are both < -180, which means the extent will be reported
2126 // incorrectly and the plypoint structure will be wrong This case is seen
2127 // first on 81004_1.KAP, (Mariannas)
2128
2129 if ((m_LonMax < -180.) && (m_LonMin < -180.)) {
2130 m_LonMin += 360.; // Normalize the extents
2131 m_LonMax += 360.;
2132
2133 Plypoint *ppp = (Plypoint *)GetCOVRTableHead(0); // Normalize the plypoints
2134 int cnPlypoint = GetCOVRTablenPoints(0);
2135
2136 for (int u = 0; u < cnPlypoint; u++) {
2137 ppp->lnp += 360.;
2138 ppp++;
2139 }
2140 }
2141
2142 return true;
2143}
2144
2145void ChartBaseBSB::SetColorScheme(ColorScheme cs, bool bApplyImmediate) {
2146 // Here we convert (subjectively) the Global ColorScheme
2147 // to an appropriate BSB_Color_Capability index.
2148
2149 switch (cs) {
2150 case GLOBAL_COLOR_SCHEME_RGB:
2151 m_mapped_color_index = COLOR_RGB_DEFAULT;
2152 break;
2153 case GLOBAL_COLOR_SCHEME_DAY:
2154 m_mapped_color_index = DAY;
2155 break;
2156 case GLOBAL_COLOR_SCHEME_DUSK:
2157 m_mapped_color_index = DUSK;
2158 break;
2159 case GLOBAL_COLOR_SCHEME_NIGHT:
2160 m_mapped_color_index = NIGHT;
2161 break;
2162 default:
2163 m_mapped_color_index = DAY;
2164 break;
2165 }
2166
2167 pPalette = GetPalettePtr(m_mapped_color_index);
2168
2169 m_global_color_scheme = cs;
2170
2171 // Force a cache dump in a simple sideways manner
2172 if (bApplyImmediate) {
2173 m_cached_scale_ppm = 1.0;
2174 }
2175
2176 // Force a new thumbnail
2177 if (pThumbData) pThumbData->pDIBThumb = NULL;
2178}
2179
2180wxBitmap *ChartBaseBSB::CreateThumbnail(int tnx, int tny, ColorScheme cs) {
2181 // Calculate the size and divisors
2182
2183 int divx = wxMax(1, Size_X / (4 * tnx));
2184 int divy = wxMax(1, Size_Y / (4 * tny));
2185
2186 int div_factor = std::min(divx, divy);
2187
2188 int des_width = Size_X / div_factor;
2189 int des_height = Size_Y / div_factor;
2190
2191 wxRect gts;
2192 gts.x = 0; // full chart
2193 gts.y = 0;
2194 gts.width = Size_X;
2195 gts.height = Size_Y;
2196
2197 int this_bpp = 24; // for wxImage
2198 // Allocate the pixel storage needed for one line of chart bits
2199 unsigned char *pLineT = (unsigned char *)malloc((Size_X + 1) * BPP / 8);
2200
2201 // Scale the data quickly
2202 unsigned char *pPixTN =
2203 (unsigned char *)malloc(des_width * des_height * this_bpp / 8);
2204
2205 int ix = 0;
2206 int iy = 0;
2207 int iyd = 0;
2208 int ixd = 0;
2209 int yoffd;
2210 unsigned char *pxs;
2211 unsigned char *pxd;
2212
2213 // Temporarily set the color scheme
2214 ColorScheme cs_tmp = m_global_color_scheme;
2215 SetColorScheme(cs, false);
2216
2217 while (iyd < des_height) {
2218 if (0 == BSBGetScanline(pLineT, iy, 0, Size_X, 1)) // get a line
2219 {
2220 free(pLineT);
2221 free(pPixTN);
2222 return NULL;
2223 }
2224
2225 yoffd = iyd * des_width * this_bpp / 8; // destination y
2226
2227 ix = 0;
2228 ixd = 0;
2229 while (ixd < des_width) {
2230 pxs = pLineT + (ix * BPP / 8);
2231 pxd = pPixTN + (yoffd + (ixd * this_bpp / 8));
2232 *pxd++ = *pxs++;
2233 *pxd++ = *pxs++;
2234 *pxd = *pxs;
2235
2236 ix += div_factor;
2237 ixd++;
2238 }
2239
2240 iy += div_factor;
2241 iyd++;
2242 }
2243
2244 free(pLineT);
2245
2246 // Reset ColorScheme
2247 SetColorScheme(cs_tmp, false);
2248
2249 wxBitmap *retBMP;
2250
2251#ifdef ocpnUSE_ocpnBitmap
2252 wxBitmap *bmx2 = new ocpnBitmap(pPixTN, des_width, des_height, -1);
2253 wxImage imgx2 = bmx2->ConvertToImage();
2254 imgx2.Rescale(des_width / 4, des_height / 4, wxIMAGE_QUALITY_HIGH);
2255 retBMP = new wxBitmap(imgx2);
2256 delete bmx2;
2257#else
2258 wxImage thumb_image(des_width, des_height, pPixTN, true);
2259 thumb_image.Rescale(des_width / 4, des_height / 4, wxIMAGE_QUALITY_HIGH);
2260 retBMP = new wxBitmap(thumb_image);
2261#endif
2262
2263 free(pPixTN);
2264
2265 return retBMP;
2266}
2267
2268//-------------------------------------------------------------------------------------------------
2269// Get the Chart thumbnail data structure
2270// Creating the thumbnail bitmap as required
2271//-------------------------------------------------------------------------------------------------
2272
2273ThumbData *ChartBaseBSB::GetThumbData(int tnx, int tny, float lat, float lon) {
2274 // Create the bitmap if needed
2275 if (!pThumbData->pDIBThumb)
2276 pThumbData->pDIBThumb = CreateThumbnail(tnx, tny, m_global_color_scheme);
2277
2278 pThumbData->Thumb_Size_X = tnx;
2279 pThumbData->Thumb_Size_Y = tny;
2280
2281 // Plot the supplied Lat/Lon on the thumbnail
2282 int divx = Size_X / tnx;
2283 int divy = Size_Y / tny;
2284
2285 int div_factor = std::min(divx, divy);
2286
2287 double pixx, pixy;
2288
2289 // Using a temporary synthetic ViewPort and source rectangle,
2290 // calculate the ships position on the thumbnail
2291 ViewPort tvp;
2292 tvp.pix_width = tnx;
2293 tvp.pix_height = tny;
2294 tvp.view_scale_ppm = GetPPM() / div_factor;
2295 wxRect trex = Rsrc;
2296 Rsrc.x = 0;
2297 Rsrc.y = 0;
2298 latlong_to_pix_vp(lat, lon, pixx, pixy, tvp);
2299 Rsrc = trex;
2300
2301 pThumbData->ShipX = pixx; // / div_factor;
2302 pThumbData->ShipY = pixy; // / div_factor;
2303
2304 return pThumbData;
2305}
2306
2307bool ChartBaseBSB::UpdateThumbData(double lat, double lon) {
2308 // Plot the supplied Lat/Lon on the thumbnail
2309 // Return TRUE if the pixel location of ownship has changed
2310
2311 int divx = Size_X / pThumbData->Thumb_Size_X;
2312 int divy = Size_Y / pThumbData->Thumb_Size_Y;
2313
2314 int div_factor = std::min(divx, divy);
2315
2316 double pixx_test, pixy_test;
2317
2318 // Using a temporary synthetic ViewPort and source rectangle,
2319 // calculate the ships position on the thumbnail
2320 ViewPort tvp;
2321 tvp.pix_width = pThumbData->Thumb_Size_X;
2322 tvp.pix_height = pThumbData->Thumb_Size_Y;
2323 tvp.view_scale_ppm = GetPPM() / div_factor;
2324 wxRect trex = Rsrc;
2325 Rsrc.x = 0;
2326 Rsrc.y = 0;
2327 latlong_to_pix_vp(lat, lon, pixx_test, pixy_test, tvp);
2328 Rsrc = trex;
2329
2330 if ((pixx_test != pThumbData->ShipX) || (pixy_test != pThumbData->ShipY)) {
2331 pThumbData->ShipX = pixx_test;
2332 pThumbData->ShipY = pixy_test;
2333 return TRUE;
2334 } else
2335 return FALSE;
2336}
2337
2338//-----------------------------------------------------------------------
2339// Pixel to Lat/Long Conversion helpers
2340//-----------------------------------------------------------------------
2341static double polytrans(double *coeff, double lon, double lat);
2342
2343int ChartBaseBSB::vp_pix_to_latlong(ViewPort &vp, double pixx, double pixy,
2344 double *plat, double *plon) {
2345 if (bHaveEmbeddedGeoref) {
2346 double raster_scale = GetPPM() / vp.view_scale_ppm;
2347
2348 double px = pixx * raster_scale + Rsrc.x;
2349 double py = pixy * raster_scale + Rsrc.y;
2350 // pix_to_latlong(px, py, plat, plon);
2351
2352 if (1) {
2353 double lon = polytrans(pwx, px, py);
2354 lon = (lon < 0) ? lon + m_cph : lon - m_cph;
2355 *plon = lon - m_lon_datum_adjust;
2356 *plat = polytrans(pwy, px, py) - m_lat_datum_adjust;
2357 }
2358
2359 return 0;
2360 } else {
2361 double slat, slon;
2362 double xp, yp;
2363
2364 if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
2365 // Use Projected Polynomial algorithm
2366
2367 double raster_scale = GetPPM() / vp.view_scale_ppm;
2368
2369 // Apply poly solution to vp center point
2370 double easting, northing;
2371 toTM(vp.clat + m_lat_datum_adjust, vp.clon + m_lon_datum_adjust,
2372 m_proj_lat, m_proj_lon, &easting, &northing);
2373 double xc = polytrans(cPoints.wpx, easting, northing);
2374 double yc = polytrans(cPoints.wpy, easting, northing);
2375
2376 // convert screen pixels to chart pixmap relative
2377 double px = xc + (pixx - (vp.pix_width / 2)) * raster_scale;
2378 double py = yc + (pixy - (vp.pix_height / 2)) * raster_scale;
2379
2380 // Apply polynomial solution to chart relative pixels to get e/n
2381 double east = polytrans(cPoints.pwx, px, py);
2382 double north = polytrans(cPoints.pwy, px, py);
2383
2384 // Apply inverse Projection to get lat/lon
2385 double lat, lon;
2386 fromTM(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2387
2388 // Datum adjustments.....
2389 //?? lon = (lon < 0) ? lon + m_cph : lon - m_cph;
2390 double slon_p = lon - m_lon_datum_adjust;
2391 double slat_p = lat - m_lat_datum_adjust;
2392
2393 // printf("%8g %8g %8g %8g %g\n", slat, slat_p, slon,
2394 // slon_p, slon - slon_p);
2395 slon = slon_p;
2396 slat = slat_p;
2397
2398 } else if (m_projection == PROJECTION_MERCATOR) {
2399 // Use Projected Polynomial algorithm
2400
2401 double raster_scale = GetPPM() / vp.view_scale_ppm;
2402
2403 // Apply poly solution to vp center point
2404 double easting, northing;
2405 toSM_ECC(vp.clat + m_lat_datum_adjust, vp.clon + m_lon_datum_adjust,
2406 m_proj_lat, m_proj_lon, &easting, &northing);
2407 double xc = polytrans(cPoints.wpx, easting, northing);
2408 double yc = polytrans(cPoints.wpy, easting, northing);
2409
2410 // convert screen pixels to chart pixmap relative
2411 double px = xc + (pixx - (vp.pix_width / 2)) * raster_scale;
2412 double py = yc + (pixy - (vp.pix_height / 2)) * raster_scale;
2413
2414 // Apply polynomial solution to chart relative pixels to get e/n
2415 double east = polytrans(cPoints.pwx, px, py);
2416 double north = polytrans(cPoints.pwy, px, py);
2417
2418 // Apply inverse Projection to get lat/lon
2419 double lat, lon;
2420 fromSM_ECC(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2421
2422 // Make Datum adjustments.....
2423 double slon_p = lon - m_lon_datum_adjust;
2424 double slat_p = lat - m_lat_datum_adjust;
2425
2426 slon = slon_p;
2427 slat = slat_p;
2428
2429 // printf("vp.clon %g xc %g px %g east %g
2430 // \n", vp.clon, xc, px, east);
2431
2432 } else if (m_projection == PROJECTION_POLYCONIC) {
2433 // Use Projected Polynomial algorithm
2434
2435 double raster_scale = GetPPM() / vp.view_scale_ppm;
2436
2437 // Apply poly solution to vp center point
2438 double easting, northing;
2439 toPOLY(vp.clat + m_lat_datum_adjust, vp.clon + m_lon_datum_adjust,
2440 m_proj_lat, m_proj_lon, &easting, &northing);
2441 double xc = polytrans(cPoints.wpx, easting, northing);
2442 double yc = polytrans(cPoints.wpy, easting, northing);
2443
2444 // convert screen pixels to chart pixmap relative
2445 double px = xc + (pixx - (vp.pix_width / 2)) * raster_scale;
2446 double py = yc + (pixy - (vp.pix_height / 2)) * raster_scale;
2447
2448 // Apply polynomial solution to chart relative pixels to get e/n
2449 double east = polytrans(cPoints.pwx, px, py);
2450 double north = polytrans(cPoints.pwy, px, py);
2451
2452 // Apply inverse Projection to get lat/lon
2453 double lat, lon;
2454 fromPOLY(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2455
2456 // Make Datum adjustments.....
2457 double slon_p = lon - m_lon_datum_adjust;
2458 double slat_p = lat - m_lat_datum_adjust;
2459
2460 slon = slon_p;
2461 slat = slat_p;
2462
2463 } else {
2464 // Use a Mercator estimator, with Eccentricity corrrection applied
2465 int dx = pixx - (vp.pix_width / 2);
2466 int dy = (vp.pix_height / 2) - pixy;
2467
2468 xp = (dx * cos(vp.skew)) - (dy * sin(vp.skew));
2469 yp = (dy * cos(vp.skew)) + (dx * sin(vp.skew));
2470
2471 double d_east = xp / vp.view_scale_ppm;
2472 double d_north = yp / vp.view_scale_ppm;
2473
2474 fromSM_ECC(d_east, d_north, vp.clat, vp.clon, &slat, &slon);
2475 }
2476
2477 *plat = slat;
2478
2479 if (slon < -180.)
2480 slon += 360.;
2481 else if (slon > 180.)
2482 slon -= 360.;
2483 *plon = slon;
2484
2485 return 0;
2486 }
2487}
2488
2489int ChartBaseBSB::latlong_to_pix_vp(double lat, double lon, double &pixx,
2490 double &pixy, ViewPort &vp) {
2491 double alat, alon;
2492
2493 if (bHaveEmbeddedGeoref) {
2494 double alat, alon;
2495
2496 alon = lon + m_lon_datum_adjust;
2497 alat = lat + m_lat_datum_adjust;
2498
2499 AdjustLongitude(alon);
2500
2501 if (1) {
2502 /* change longitude phase (CPH) */
2503 double lonp = (alon < 0) ? alon + m_cph : alon - m_cph;
2504 double xd = polytrans(wpx, lonp, alat);
2505 double yd = polytrans(wpy, lonp, alat);
2506
2507 double raster_scale = GetPPM() / vp.view_scale_ppm;
2508
2509 pixx = (xd - Rsrc.x) / raster_scale;
2510 pixy = (yd - Rsrc.y) / raster_scale;
2511
2512 return 0;
2513 }
2514 } else {
2515 double easting, northing;
2516 double xlon = lon;
2517
2518 // Make sure lon and lon0 are same phase
2519 /*
2520 if((xlon * vp.clon) < 0.)
2521 {
2522 if(xlon < 0.)
2523 xlon += 360.;
2524 else
2525 xlon -= 360.;
2526 }
2527
2528 if(fabs(xlon - vp.clon) > 180.)
2529 {
2530 if(xlon > vp.clon)
2531 xlon -= 360.;
2532 else
2533 xlon += 360.;
2534 }
2535 */
2536
2537 if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
2538 // Use Projected Polynomial algorithm
2539
2540 alon = lon + m_lon_datum_adjust;
2541 alat = lat + m_lat_datum_adjust;
2542
2543 // Get e/n from TM Projection
2544 toTM(alat, alon, m_proj_lat, m_proj_lon, &easting, &northing);
2545
2546 // Apply poly solution to target point
2547 double xd = polytrans(cPoints.wpx, easting, northing);
2548 double yd = polytrans(cPoints.wpy, easting, northing);
2549
2550 // Apply poly solution to vp center point
2551 toTM(vp.clat + m_lat_datum_adjust, vp.clon + m_lon_datum_adjust,
2552 m_proj_lat, m_proj_lon, &easting, &northing);
2553 double xc = polytrans(cPoints.wpx, easting, northing);
2554 double yc = polytrans(cPoints.wpy, easting, northing);
2555
2556 // Calculate target point relative to vp center
2557 double raster_scale = GetPPM() / vp.view_scale_ppm;
2558
2559 double xs = xc - vp.pix_width * raster_scale / 2;
2560 double ys = yc - vp.pix_height * raster_scale / 2;
2561
2562 pixx = (xd - xs) / raster_scale;
2563 pixy = (yd - ys) / raster_scale;
2564
2565 } else if (m_projection == PROJECTION_MERCATOR) {
2566 // Use Projected Polynomial algorithm
2567
2568 alon = lon + m_lon_datum_adjust;
2569 alat = lat + m_lat_datum_adjust;
2570
2571 // Get e/n from Projection
2572 xlon = alon;
2573 AdjustLongitude(xlon);
2574 toSM_ECC(alat, xlon, m_proj_lat, m_proj_lon, &easting, &northing);
2575
2576 // Apply poly solution to target point
2577 double xd = polytrans(cPoints.wpx, easting, northing);
2578 double yd = polytrans(cPoints.wpy, easting, northing);
2579
2580 // Apply poly solution to vp center point
2581 double xlonc = vp.clon;
2582 AdjustLongitude(xlonc);
2583
2584 toSM_ECC(vp.clat + m_lat_datum_adjust, xlonc + m_lon_datum_adjust,
2585 m_proj_lat, m_proj_lon, &easting, &northing);
2586 double xc = polytrans(cPoints.wpx, easting, northing);
2587 double yc = polytrans(cPoints.wpy, easting, northing);
2588
2589 // Calculate target point relative to vp center
2590 double raster_scale = GetPPM() / vp.view_scale_ppm;
2591
2592 double xs = xc - vp.pix_width * raster_scale / 2;
2593 double ys = yc - vp.pix_height * raster_scale / 2;
2594
2595 pixx = (xd - xs) / raster_scale;
2596 pixy = (yd - ys) / raster_scale;
2597
2598 } else if (m_projection == PROJECTION_POLYCONIC) {
2599 // Use Projected Polynomial algorithm
2600
2601 alon = lon + m_lon_datum_adjust;
2602 alat = lat + m_lat_datum_adjust;
2603
2604 // Get e/n from Projection
2605 xlon = AdjustLongitude(alon);
2606 toPOLY(alat, xlon, m_proj_lat, m_proj_lon, &easting, &northing);
2607
2608 // Apply poly solution to target point
2609 double xd = polytrans(cPoints.wpx, easting, northing);
2610 double yd = polytrans(cPoints.wpy, easting, northing);
2611
2612 // Apply poly solution to vp center point
2613 double xlonc = AdjustLongitude(vp.clon);
2614
2615 toPOLY(vp.clat + m_lat_datum_adjust, xlonc + m_lon_datum_adjust,
2616 m_proj_lat, m_proj_lon, &easting, &northing);
2617 double xc = polytrans(cPoints.wpx, easting, northing);
2618 double yc = polytrans(cPoints.wpy, easting, northing);
2619
2620 // Calculate target point relative to vp center
2621 double raster_scale = GetPPM() / vp.view_scale_ppm;
2622
2623 double xs = xc - vp.pix_width * raster_scale / 2;
2624 double ys = yc - vp.pix_height * raster_scale / 2;
2625
2626 pixx = (xd - xs) / raster_scale;
2627 pixy = (yd - ys) / raster_scale;
2628
2629 } else {
2630 toSM_ECC(lat, xlon, vp.clat, vp.clon, &easting, &northing);
2631
2632 double epix = easting * vp.view_scale_ppm;
2633 double npix = northing * vp.view_scale_ppm;
2634
2635 double dx = epix * cos(vp.skew) + npix * sin(vp.skew);
2636 double dy = npix * cos(vp.skew) - epix * sin(vp.skew);
2637
2638 pixx = ((double)vp.pix_width / 2) + dx;
2639 pixy = ((double)vp.pix_height / 2) - dy;
2640 }
2641 return 0;
2642 }
2643
2644 return 1;
2645}
2646
2647void ChartBaseBSB::latlong_to_chartpix(double lat, double lon, double &pixx,
2648 double &pixy) {
2649 double alat, alon;
2650 pixx = 0.0;
2651 pixy = 0.0;
2652
2653 if (bHaveEmbeddedGeoref) {
2654 double alat, alon;
2655
2656 alon = lon + m_lon_datum_adjust;
2657 alat = lat + m_lat_datum_adjust;
2658
2659 alon = AdjustLongitude(alon);
2660
2661 /* change longitude phase (CPH) */
2662 double lonp = (alon < 0) ? alon + m_cph : alon - m_cph;
2663 pixx = polytrans(wpx, lonp, alat);
2664 pixy = polytrans(wpy, lonp, alat);
2665 } else {
2666 double easting, northing;
2667 double xlon = lon;
2668
2669 if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
2670 // Use Projected Polynomial algorithm
2671
2672 alon = lon + m_lon_datum_adjust;
2673 alat = lat + m_lat_datum_adjust;
2674
2675 // Get e/n from TM Projection
2676 toTM(alat, alon, m_proj_lat, m_proj_lon, &easting, &northing);
2677
2678 // Apply poly solution to target point
2679 pixx = polytrans(cPoints.wpx, easting, northing);
2680 pixy = polytrans(cPoints.wpy, easting, northing);
2681
2682 } else if (m_projection == PROJECTION_MERCATOR) {
2683 // Use Projected Polynomial algorithm
2684
2685 alon = lon + m_lon_datum_adjust;
2686 alat = lat + m_lat_datum_adjust;
2687
2688 // Get e/n from Projection
2689 xlon = AdjustLongitude(alon);
2690
2691 toSM_ECC(alat, xlon, m_proj_lat, m_proj_lon, &easting, &northing);
2692
2693 // Apply poly solution to target point
2694 pixx = polytrans(cPoints.wpx, easting, northing);
2695 pixy = polytrans(cPoints.wpy, easting, northing);
2696
2697 } else if (m_projection == PROJECTION_POLYCONIC) {
2698 // Use Projected Polynomial algorithm
2699
2700 alon = lon + m_lon_datum_adjust;
2701 alat = lat + m_lat_datum_adjust;
2702
2703 // Get e/n from Projection
2704 xlon = AdjustLongitude(alon);
2705 toPOLY(alat, xlon, m_proj_lat, m_proj_lon, &easting, &northing);
2706
2707 // Apply poly solution to target point
2708 pixx = polytrans(cPoints.wpx, easting, northing);
2709 pixy = polytrans(cPoints.wpy, easting, northing);
2710 }
2711 }
2712}
2713
2714void ChartBaseBSB::chartpix_to_latlong(double pixx, double pixy, double *plat,
2715 double *plon) {
2716 if (bHaveEmbeddedGeoref) {
2717 double lon = polytrans(pwx, pixx, pixy);
2718 lon = (lon < 0) ? lon + m_cph : lon - m_cph;
2719 *plon = lon - m_lon_datum_adjust;
2720 *plat = polytrans(pwy, pixx, pixy) - m_lat_datum_adjust;
2721 } else {
2722 double slat, slon;
2723 if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
2724 // Use Projected Polynomial algorithm
2725
2726 // Apply polynomial solution to chart relative pixels to get e/n
2727 double east = polytrans(cPoints.pwx, pixx, pixy);
2728 double north = polytrans(cPoints.pwy, pixx, pixy);
2729
2730 // Apply inverse Projection to get lat/lon
2731 double lat, lon;
2732 fromTM(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2733
2734 // Datum adjustments.....
2735 //?? lon = (lon < 0) ? lon + m_cph : lon - m_cph;
2736 slon = lon - m_lon_datum_adjust;
2737 slat = lat - m_lat_datum_adjust;
2738
2739 } else if (m_projection == PROJECTION_MERCATOR) {
2740 // Use Projected Polynomial algorithm
2741 // Apply polynomial solution to chart relative pixels to get e/n
2742 double east = polytrans(cPoints.pwx, pixx, pixy);
2743 double north = polytrans(cPoints.pwy, pixx, pixy);
2744
2745 // Apply inverse Projection to get lat/lon
2746 double lat, lon;
2747 fromSM_ECC(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2748
2749 // Make Datum adjustments.....
2750 slon = lon - m_lon_datum_adjust;
2751 slat = lat - m_lat_datum_adjust;
2752 } else if (m_projection == PROJECTION_POLYCONIC) {
2753 // Use Projected Polynomial algorithm
2754 // Apply polynomial solution to chart relative pixels to get e/n
2755 double east = polytrans(cPoints.pwx, pixx, pixy);
2756 double north = polytrans(cPoints.pwy, pixx, pixy);
2757
2758 // Apply inverse Projection to get lat/lon
2759 double lat, lon;
2760 fromPOLY(east, north, m_proj_lat, m_proj_lon, &lat, &lon);
2761
2762 // Make Datum adjustments.....
2763 slon = lon - m_lon_datum_adjust;
2764 slat = lat - m_lat_datum_adjust;
2765
2766 } else {
2767 slon = 0.;
2768 slat = 0.;
2769 }
2770
2771 *plat = slat;
2772
2773 if (slon < -180.)
2774 slon += 360.;
2775 else if (slon > 180.)
2776 slon -= 360.;
2777 *plon = slon;
2778 }
2779}
2780
2781void ChartBaseBSB::ComputeSourceRectangle(const ViewPort &vp,
2782 wxRect *pSourceRect) {
2783 m_raster_scale_factor = GetRasterScaleFactor(vp);
2784 double xd, yd;
2785 latlong_to_chartpix(vp.clat, vp.clon, xd, yd);
2786
2787 wxRealPoint pos, size;
2788
2789 pos.x = xd - (vp.pix_width * m_raster_scale_factor / 2);
2790 pos.y = yd - (vp.pix_height * m_raster_scale_factor / 2);
2791
2792 size.x = vp.pix_width * m_raster_scale_factor;
2793 size.y = vp.pix_height * m_raster_scale_factor;
2794
2795 *pSourceRect =
2796 wxRect(wxRound(pos.x), wxRound(pos.y), wxRound(size.x), wxRound(size.y));
2797}
2798
2799double ChartBaseBSB::GetRasterScaleFactor(const ViewPort &vp) {
2800 // This funny contortion is necessary to allow scale factors < 1, i.e.
2801 // overzoom
2802 return (wxRound(100000 * GetPPM() / vp.view_scale_ppm)) / 100000.;
2803}
2804
2805void ChartBaseBSB::SetVPRasterParms(const ViewPort &vpt) {
2806 // Calculate the potential datum offset parameters for this viewport, if
2807 // not WGS84
2808
2809 if (m_datum_index == DATUM_INDEX_WGS84 ||
2810 m_datum_index == DATUM_INDEX_UNKNOWN) {
2811 m_lon_datum_adjust = (-m_dtm_lon) / 3600.;
2812 m_lat_datum_adjust = (-m_dtm_lat) / 3600.;
2813 }
2814
2815 else {
2816 double to_lat, to_lon;
2817 MolodenskyTransform(vpt.clat, vpt.clon, &to_lat, &to_lon, m_datum_index,
2818 DATUM_INDEX_WGS84);
2819 m_lon_datum_adjust = -(to_lon - vpt.clon);
2820 m_lat_datum_adjust = -(to_lat - vpt.clat);
2821 if (m_b_apply_dtm) {
2822 m_lon_datum_adjust -= m_dtm_lon / 3600.;
2823 m_lat_datum_adjust -= m_dtm_lat / 3600.;
2824 }
2825 }
2826
2827 ComputeSourceRectangle(vpt, &Rsrc);
2828
2829 if (vpt.IsValid()) m_vp_render_last = vpt;
2830}
2831
2832bool ChartBaseBSB::AdjustVP(ViewPort &vp_last, ViewPort &vp_proposed) {
2833 bool bInside = G_FloatPtInPolygon((MyFlPoint *)GetCOVRTableHead(0),
2834 GetCOVRTablenPoints(0), vp_proposed.clon,
2835 vp_proposed.clat);
2836 if (!bInside) return false;
2837
2838 int ret_val = 0;
2839 double binary_scale_factor = GetPPM() / vp_proposed.view_scale_ppm;
2840
2841 if (vp_last.IsValid()) {
2842 // We only need to adjust the VP if the cache is valid and potentially
2843 // usable, i.e. the scale factor is integer... The objective here is to
2844 // ensure that the VP center falls on an exact pixel boundary within the
2845 // cache
2846
2847 if (cached_image_ok && (binary_scale_factor > 1.0) &&
2848 (fabs(binary_scale_factor - wxRound(binary_scale_factor)) < 1e-5)) {
2849 if (m_b_cdebug)
2850 printf(" Possible Adjust VP for integer scale: %g\n",
2851 binary_scale_factor);
2852
2853 wxRect rprop;
2854 ComputeSourceRectangle(vp_proposed, &rprop);
2855
2856 double pixx, pixy;
2857 double lon_adj, lat_adj;
2858 latlong_to_pix_vp(vp_proposed.clat, vp_proposed.clon, pixx, pixy,
2859 vp_proposed);
2860 vp_pix_to_latlong(vp_proposed, pixx, pixy, &lat_adj, &lon_adj);
2861
2862 vp_proposed.clat = lat_adj;
2863 vp_proposed.clon = lon_adj;
2864 ret_val = 1;
2865 }
2866 }
2867
2868 return (ret_val > 0);
2869}
2870
2871bool ChartBaseBSB::IsRenderCacheable(wxRect &source, wxRect &dest) {
2872 double scale_x = (double)source.width / (double)dest.width;
2873
2874 if (scale_x <= 1.0) // overzoom
2875 {
2876 // if(m_b_cdebug)printf(" MISS<<<>>>GVUC: Overzoom\n");
2877 return false;
2878 }
2879
2880 // Using the cache only works for pure binary scale factors......
2881 if ((fabs(scale_x - wxRound(scale_x))) > .0001) {
2882 // if(m_b_cdebug)printf(" MISS<<<>>>GVUC: Not digital scale
2883 // test 1\n");
2884 return false;
2885 }
2886
2887 // Scale must be exactly digital...
2888 if ((int)(source.width / dest.width) != (int)wxRound(scale_x)) {
2889 // if(m_b_cdebug)printf(" MISS<<<>>>GVUC: Not digital scale
2890 // test 2\n");
2891 return false;
2892 }
2893
2894 return true;
2895}
2896
2897void ChartBaseBSB::GetValidCanvasRegion(const ViewPort &VPoint,
2898 OCPNRegion *pValidRegion) {
2899 SetVPRasterParms(VPoint);
2900
2901 double raster_scale = VPoint.view_scale_ppm / GetPPM();
2902
2903 int rxl, rxr;
2904 int ryb, ryt;
2905
2906 rxl = wxMax(-Rsrc.x * raster_scale, VPoint.rv_rect.x);
2907 rxr = wxMin((Size_X - Rsrc.x) * raster_scale,
2908 VPoint.rv_rect.width + VPoint.rv_rect.x);
2909
2910 ryt = wxMax(-Rsrc.y * raster_scale, VPoint.rv_rect.y);
2911 ryb = wxMin((Size_Y - Rsrc.y) * raster_scale,
2912 VPoint.rv_rect.height + VPoint.rv_rect.y);
2913
2914 pValidRegion->Clear();
2915 pValidRegion->Union(rxl, ryt, rxr - rxl, ryb - ryt);
2916}
2917
2918LLRegion ChartBaseBSB::GetValidRegion() {
2919 // should we cache this?
2920 double ll[8];
2921 chartpix_to_latlong(0, 0, ll + 0, ll + 1);
2922 chartpix_to_latlong(0, Size_Y, ll + 2, ll + 3);
2923 chartpix_to_latlong(Size_X, Size_Y, ll + 4, ll + 5);
2924 chartpix_to_latlong(Size_X, 0, ll + 6, ll + 7);
2925
2926 // for now don't allow raster charts to cross both 0 meridian and IDL
2927 // (complicated to deal with)
2928 for (int i = 1; i < 6; i += 2)
2929 if (fabs(ll[i] - ll[i + 2]) > 180) {
2930 // we detect crossing idl here, make all longitudes positive
2931 for (int i = 1; i < 8; i += 2)
2932 if (ll[i] < 0) ll[i] += 360;
2933 break;
2934 }
2935
2936 return LLRegion(4, ll);
2937}
2938
2939bool ChartBaseBSB::GetViewUsingCache(wxRect &source, wxRect &dest,
2940 const OCPNRegion &Region,
2941 ScaleTypeEnum scale_type) {
2942 wxRect s1;
2943 ScaleTypeEnum scale_type_corrected;
2944
2945 if (m_b_cdebug) printf(" source: %d %d\n", source.x, source.y);
2946 if (m_b_cdebug) printf(" cache: %d %d\n", cache_rect.x, cache_rect.y);
2947
2948 // Anything to do?
2949 if ((source == cache_rect) /*&& (cache_scale_method == scale_type)*/ &&
2950 (cached_image_ok)) {
2951 if (m_b_cdebug) printf(" GVUC: Cache is good, nothing to do\n");
2952 return false;
2953 }
2954
2955 double scale_x = (double)source.width / (double)dest.width;
2956
2957 if (m_b_cdebug) printf("GVUC: scale_x: %g\n", scale_x);
2958
2959 // Enforce a limit on bilinear scaling, for performance reasons
2960 scale_type_corrected = scale_type; // RENDER_LODEF; //scale_type;
2961 if (scale_x > m_bilinear_limit) scale_type_corrected = RENDER_LODEF;
2962
2963 {
2964 // if(b_cdebug)printf(" MISS<<<>>>GVUC: Intentional out\n");
2965 // return GetView( source, dest, scale_type_corrected );
2966 }
2967
2968 // Using the cache only works for pure binary scale factors......
2969 if ((fabs(scale_x - wxRound(scale_x))) > .0001) {
2970 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: Not digital scale test 1\n");
2971 return GetView(source, dest, scale_type_corrected);
2972 }
2973
2974 // scale_type_corrected = RENDER_LODEF;
2975
2976 if (!cached_image_ok) {
2977 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: Cache NOk\n");
2978 return GetView(source, dest, scale_type_corrected);
2979 }
2980
2981 if (scale_x < 1.0) // overzoom
2982 {
2983 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: Overzoom\n");
2984 return GetView(source, dest, scale_type_corrected);
2985 }
2986
2987 // Scale must be exactly digital...
2988 if ((int)(source.width / dest.width) != (int)wxRound(scale_x)) {
2989 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: Not digital scale test 2\n");
2990 return GetView(source, dest, scale_type_corrected);
2991 }
2992
2993 // Calculate the digital scale, e.g. 1,2,4,8,,,
2994 int cs1d = source.width / dest.width;
2995 if (abs(source.x - cache_rect.x) % cs1d) {
2996 if (m_b_cdebug)
2997 printf(" source.x: %d cache_rect.x: %d cs1d: %d\n", source.x,
2998 cache_rect.x, cs1d);
2999 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: x mismatch\n");
3000 return GetView(source, dest, scale_type_corrected);
3001 }
3002 if (abs(source.y - cache_rect.y) % cs1d) {
3003 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: y mismatch\n");
3004 return GetView(source, dest, scale_type_corrected);
3005 }
3006
3007 if (pPixCache && ((pPixCache->GetWidth() != dest.width) ||
3008 (pPixCache->GetHeight() != dest.height))) {
3009 if (m_b_cdebug) printf(" MISS<<<>>>GVUC: dest size mismatch\n");
3010 return GetView(source, dest, scale_type_corrected);
3011 }
3012
3013 int stride_rows =
3014 (source.y + source.height) - (cache_rect.y + cache_rect.height);
3015 int scaled_stride_rows = (int)(stride_rows / scale_x);
3016
3017 if (abs(stride_rows) >= source.height) // Pan more than one screen
3018 return GetView(source, dest, scale_type_corrected);
3019
3020 int stride_pixels =
3021 (source.x + source.width) - (cache_rect.x + cache_rect.width);
3022 int scaled_stride_pixels = (int)(stride_pixels / scale_x);
3023
3024 if (abs(stride_pixels) >= source.width) // Pan more than one screen
3025 return GetView(source, dest, scale_type_corrected);
3026
3027 if (m_b_cdebug) printf(" GVUC Using raster data cache\n");
3028
3029 ScaleTypeEnum pan_scale_type_x = scale_type_corrected;
3030 ScaleTypeEnum pan_scale_type_y = scale_type_corrected;
3031
3032 // "Blit" the valid pixels out of the way
3033 if (pPixCache) {
3034 int height = pPixCache->GetHeight();
3035 int width = pPixCache->GetWidth();
3036 int buffer_stride_bytes = pPixCache->GetLinePitch();
3037
3038 unsigned char *ps;
3039 unsigned char *pd;
3040
3041 if (stride_rows > 0) // pan down
3042 {
3043 ps = pPixCache->GetpData() +
3044 (abs(scaled_stride_rows) * buffer_stride_bytes);
3045 if (stride_pixels > 0) ps += scaled_stride_pixels * BPP / 8;
3046
3047 pd = pPixCache->GetpData();
3048 if (stride_pixels <= 0) pd += abs(scaled_stride_pixels) * BPP / 8;
3049
3050 for (int iy = 0; iy < (height - abs(scaled_stride_rows)); iy++) {
3051 memmove(pd, ps, (width - abs(scaled_stride_pixels)) * BPP / 8);
3052 ps += buffer_stride_bytes;
3053 pd += buffer_stride_bytes;
3054 }
3055
3056 } else {
3057 ps = pPixCache->GetpData() +
3058 ((height - abs(scaled_stride_rows) - 1) * buffer_stride_bytes);
3059 if (stride_pixels > 0) // make a hole on right
3060 ps += scaled_stride_pixels * BPP / 8;
3061
3062 pd = pPixCache->GetpData() + ((height - 1) * buffer_stride_bytes);
3063 if (stride_pixels <= 0) // make a hole on the left
3064 pd += abs(scaled_stride_pixels) * BPP / 8;
3065
3066 for (int iy = 0; iy < (height - abs(scaled_stride_rows)); iy++) {
3067 memmove(pd, ps, (width - abs(scaled_stride_pixels)) * BPP / 8);
3068 ps -= buffer_stride_bytes;
3069 pd -= buffer_stride_bytes;
3070 }
3071 }
3072
3073 // Y Pan
3074 if (source.y != cache_rect.y) {
3075 wxRect sub_dest = dest;
3076 sub_dest.height = abs(scaled_stride_rows);
3077
3078 if (stride_rows > 0) // pan down
3079 {
3080 sub_dest.y = height - scaled_stride_rows;
3081
3082 } else {
3083 sub_dest.y = 0;
3084 }
3085
3086 // Get the new bits needed
3087
3088 // A little optimization...
3089 // No sense in fetching bits that are not part of the ultimate render
3090 // region
3091 wxRegionContain rc = Region.Contains(sub_dest);
3092 if ((wxPartRegion == rc) || (wxInRegion == rc)) {
3093 GetAndScaleData(pPixCache->GetpData(), pPixCache->GetLength(), source,
3094 source.width, sub_dest, width, cs1d, pan_scale_type_y);
3095 }
3096 pPixCache->Update();
3097
3098 // Update the cached parameters, Y only
3099
3100 cache_rect.y = source.y;
3101 // cache_rect = source;
3102 cache_rect_scaled = dest;
3103 cached_image_ok = 1;
3104
3105 } // Y Pan
3106
3107 // X Pan
3108 if (source.x != cache_rect.x) {
3109 wxRect sub_dest = dest;
3110 sub_dest.width = abs(scaled_stride_pixels);
3111
3112 if (stride_pixels > 0) // pan right
3113 {
3114 sub_dest.x = width - scaled_stride_pixels;
3115 } else // pan left
3116 {
3117 sub_dest.x = 0;
3118 }
3119
3120 // Get the new bits needed
3121
3122 // A little optimization...
3123 // No sense in fetching bits that are not part of the ultimate render
3124 // region
3125 wxRegionContain rc = Region.Contains(sub_dest);
3126 if ((wxPartRegion == rc) || (wxInRegion == rc)) {
3127 GetAndScaleData(pPixCache->GetpData(), pPixCache->GetLength(), source,
3128 source.width, sub_dest, width, cs1d, pan_scale_type_x);
3129 }
3130
3131 pPixCache->Update();
3132
3133 // Update the cached parameters
3134 cache_rect = source;
3135 cache_rect_scaled = dest;
3136 cached_image_ok = 1;
3137
3138 } // X pan
3139
3140 return true;
3141 }
3142 return false;
3143}
3144
3145int s_dc;
3146
3147bool ChartBaseBSB::RenderViewOnDC(wxMemoryDC &dc, const ViewPort &VPoint) {
3148 SetVPRasterParms(VPoint);
3149
3150 OCPNRegion rgn(0, 0, VPoint.pix_width, VPoint.pix_height);
3151
3152 bool bsame_region = (rgn == m_last_region); // only want to do this once
3153
3154 if (!bsame_region) cached_image_ok = false;
3155
3156 m_last_region = rgn;
3157
3158 return RenderRegionViewOnDC(dc, VPoint, rgn);
3159}
3160
3161bool ChartBaseBSB::RenderRegionViewOnGL(const wxGLContext &glc,
3162 const ViewPort &VPoint,
3163 const OCPNRegion &RectRegion,
3164 const LLRegion &Region) {
3165 return true;
3166}
3167
3168bool ChartBaseBSB::RenderRegionViewOnDC(wxMemoryDC &dc, const ViewPort &VPoint,
3169 const OCPNRegion &Region) {
3170 SetVPRasterParms(VPoint);
3171
3172 wxRect dest(0, 0, VPoint.pix_width, VPoint.pix_height);
3173 // double factor = ((double)Rsrc.width)/((double)dest.width);
3174 double factor = GetRasterScaleFactor(VPoint);
3175 if (m_b_cdebug) {
3176 printf("%d RenderRegion ScaleType: %d factor: %g\n", s_dc++,
3177 RENDER_HIDEF, factor);
3178 printf("Rect list:\n");
3179 OCPNRegionIterator upd(Region); // get the requested rect list
3180 while (upd.HaveRects()) {
3181 wxRect rect = upd.GetRect();
3182 printf(" %d %d %d %d\n", rect.x, rect.y, rect.width, rect.height);
3183 upd.NextRect();
3184 }
3185 }
3186
3187 // Invalidate the cache if the scale has changed or the viewport size has
3188 // changed....
3189 if ((fabs(m_cached_scale_ppm - VPoint.view_scale_ppm) > 1e-9) ||
3190 (m_last_vprect != dest)) {
3191 cached_image_ok = false;
3192 m_vp_render_last.Invalidate();
3193 }
3194 /*
3195 if(pPixCache)
3196 {
3197 if((pPixCache->GetWidth() != dest.width) ||
3198 (pPixCache->GetHeight() != dest.height))
3199 {
3200 delete pPixCache;
3201 pPixCache = new PixelCache(dest.width, dest.height, BPP);
3202 }
3203 }
3204 else
3205 pPixCache = new PixelCache(dest.width, dest.height, BPP);
3206 */
3207
3208 m_cached_scale_ppm = VPoint.view_scale_ppm;
3209 m_last_vprect = dest;
3210
3211 if (cached_image_ok) {
3212 // Anything to do?
3213 bool bsame_region = (Region == m_last_region); // only want to do this once
3214
3215 if ((bsame_region) && (Rsrc == cache_rect)) {
3216 pPixCache->SelectIntoDC(dc);
3217 if (m_b_cdebug) printf(" Using Current PixelCache\n");
3218 return false;
3219 }
3220 }
3221
3222 m_last_region = Region;
3223
3224 // Analyze the region requested
3225 // When rendering complex regions, (more than say 4 rectangles)
3226 // .OR. small proportions, then rectangle rendering may be faster
3227 // Also true if the scale is less than near unity, or overzoom.
3228 // This will be the case for backgrounds of the quilt.
3229
3230 /* Update for Version 2.4.0
3231 This logic seems flawed, at least for quilts which contain charts having
3232 non-rectangular coverage areas. These quilt regions decompose to ...LOTS... of
3233 rectangles, most of which are 1 pixel in height. This is very slow, due to the
3234 overhead of GetAndScaleData(). However, remember that overzoom never uses the
3235 cache, nor does non-binary scale factors.. So, we check to see if this is a
3236 cacheable render, and that the number of rectangles is "reasonable"
3237 */
3238
3239 // Get the region rectangle count
3240
3241 int n_rect = 0;
3242 OCPNRegionIterator upd(Region); // get the requested rect list
3243 while (upd.HaveRects()) {
3244 n_rect++;
3245 upd.NextRect();
3246 }
3247
3248 if ((!IsRenderCacheable(Rsrc, dest) && (n_rect > 4) && (n_rect < 20)) ||
3249 (factor < 1)) {
3250 if (m_b_cdebug)
3251 printf(" RenderRegion by rect iterator n_rect: %d\n", n_rect);
3252
3253 // Verify that the persistent pixel cache is at least as large as the
3254 // largest rectangle in the region
3255 wxRect dest_check_rect = dest;
3256 OCPNRegionIterator upd_check(Region); // get the requested rect list
3257 while (upd_check.HaveRects()) {
3258 wxRect rect = upd_check.GetRect();
3259 dest_check_rect.Union(rect);
3260 upd_check.NextRect();
3261 }
3262
3263 if (pPixCache) {
3264 if ((pPixCache->GetWidth() != dest_check_rect.width) ||
3265 (pPixCache->GetHeight() != dest_check_rect.height)) {
3266 delete pPixCache;
3267 pPixCache =
3268 new PixelCache(dest_check_rect.width, dest_check_rect.height, BPP);
3269 }
3270 } else
3271 pPixCache =
3272 new PixelCache(dest_check_rect.width, dest_check_rect.height, BPP);
3273
3274 ScaleTypeEnum ren_type = RENDER_LODEF;
3275
3276 // Decompose the region into rectangles, and fetch them into the target
3277 // dc
3278 OCPNRegionIterator upd(Region); // get the requested rect list
3279 int ir = 0;
3280 while (upd.HaveRects()) {
3281 wxRect rect = upd.GetRect();
3282
3283 // Floating point math can lead to negative rectangle origin.
3284 // If this happens, we arbitrarily shift the rectangle to be positive
3285 // semidefinite. This will cause at most a 1 pixlel error onscreen.
3286 if (rect.y < 0) rect.Offset(0, -rect.y);
3287 if (rect.x < 0) rect.Offset(-rect.x, 0);
3288
3289 GetAndScaleData(pPixCache->GetpData(), pPixCache->GetLength(), Rsrc,
3290 Rsrc.width, rect, pPixCache->GetWidth(), factor,
3291 ren_type);
3292
3293 ir++;
3294 upd.NextRect();
3295 ;
3296 }
3297
3298 pPixCache->Update();
3299
3300 // Update cache parameters
3301 cache_rect = Rsrc;
3302 cache_scale_method = ren_type;
3303 cached_image_ok = false; // Never cache this type of render
3304
3305 // Select the data into the dc
3306 pPixCache->SelectIntoDC(dc);
3307
3308 return true;
3309 }
3310
3311 // Default is to try using the cache
3312
3313 if (pPixCache) {
3314 if ((pPixCache->GetWidth() != dest.width) ||
3315 (pPixCache->GetHeight() != dest.height)) {
3316 delete pPixCache;
3317 pPixCache = new PixelCache(dest.width, dest.height, BPP);
3318 }
3319 } else
3320 pPixCache = new PixelCache(dest.width, dest.height, BPP);
3321
3322 if (m_b_cdebug) printf(" Render Region By GVUC\n");
3323
3324 // A performance enhancement.....
3325 ScaleTypeEnum scale_type_zoom = RENDER_HIDEF;
3326 double binary_scale_factor = VPoint.view_scale_ppm / GetPPM();
3327 if (binary_scale_factor < .20) scale_type_zoom = RENDER_LODEF;
3328
3329 bool bnewview = GetViewUsingCache(Rsrc, dest, Region, scale_type_zoom);
3330
3331 // Select the data into the dc
3332 pPixCache->SelectIntoDC(dc);
3333
3334 return bnewview;
3335}
3336
3337wxImage *ChartBaseBSB::GetImage() {
3338 int img_size_x = ((Size_X >> 2) * 4) + 4;
3339 wxImage *img = new wxImage(img_size_x, Size_Y, false);
3340
3341 unsigned char *ppnx = img->GetData();
3342
3343 for (int i = 0; i < Size_Y; i++) {
3344 wxRect source_rect(0, i, Size_X, 1);
3345 wxRect dest_rect(0, 0, Size_X, 1);
3346
3347 GetAndScaleData(img->GetData(), img_size_x * Size_Y * 3, source_rect,
3348 Size_X, dest_rect, Size_X, 1.0, RENDER_HIDEF);
3349
3350 ppnx += img_size_x * 3;
3351 }
3352
3353 return img;
3354}
3355
3356bool ChartBaseBSB::GetView(wxRect &source, wxRect &dest,
3357 ScaleTypeEnum scale_type) {
3358 assert(pPixCache != 0);
3359 // PixelCache *pPixCacheTemp = new PixelCache(dest.width, dest.height,
3360 // BPP);
3361
3362 // Get and Rescale the data directly into the temporary PixelCache data
3363 // buffer
3364 double factor = ((double)source.width) / ((double)dest.width);
3365
3366 /*
3367 if(!GetAndScaleData(&ppnx, source, source.width, dest, dest.width,
3368 factor, scale_type))
3369 {
3370 delete pPixCacheTemp; // Some error, retain
3371 old cache return false;
3372 }
3373 else
3374 {
3375 delete pPixCache; // new cache is OK
3376 pPixCache = pPixCacheTemp;
3377 }
3378 */
3379 GetAndScaleData(pPixCache->GetpData(), pPixCache->GetLength(), source,
3380 source.width, dest, dest.width, factor, scale_type);
3381 pPixCache->Update();
3382
3383 // Update cache parameters
3384
3385 cache_rect = source;
3386 cache_rect_scaled = dest;
3387 cache_scale_method = scale_type;
3388
3389 cached_image_ok = 1;
3390
3391 return TRUE;
3392}
3393
3394bool ChartBaseBSB::GetAndScaleData(unsigned char *ppn, size_t data_size,
3395 wxRect &source, int source_stride,
3396 wxRect &dest, int dest_stride,
3397 double scale_factor,
3398 ScaleTypeEnum scale_type) {
3399 unsigned char *s_data = NULL;
3400
3401 double factor = scale_factor;
3402 int Factor = (int)factor;
3403
3404 int target_width = (int)wxRound((double)source.width / factor);
3405 int target_height = (int)wxRound((double)source.height / factor);
3406
3407 int dest_line_length = dest_stride * BPP / 8;
3408
3409 // On MSW, if using DibSections, each scan line starts on a DWORD boundary.
3410 // The DibSection has been allocated to conform with this requirement.
3411#ifdef __PIX_CACHE_DIBSECTION__
3412 dest_line_length = (((dest_stride * 24) + 31) & ~31) >> 3;
3413#endif
3414
3415 if ((target_height == 0) || (target_width == 0)) return false;
3416
3417 // `volatile` may be needed with regard to `sigsetjmp()` below;
3418 // at least it prevents compiler warning about clobbering the variable.
3419 unsigned char *volatile target_data = ppn;
3420 unsigned char *data = ppn;
3421
3422 if (factor > 1) // downsampling
3423 {
3424 if (scale_type == RENDER_HIDEF) {
3425 // Allocate a working buffer based on scale factor
3426 int blur_factor = wxMax(2, Factor);
3427 int wb_size = (source.width) * (blur_factor * 2) * BPP / 8;
3428 s_data = (unsigned char *)malloc(wb_size); // work buffer
3429 unsigned char *pixel;
3430 int y_offset;
3431
3432 for (int y = dest.y; y < (dest.y + dest.height); y++) {
3433 // Read "blur_factor" lines
3434
3435 wxRect s1;
3436 s1.x = source.x;
3437 s1.y = source.y + (int)(y * factor);
3438 s1.width = source.width;
3439 s1.height = blur_factor;
3440 GetChartBits(s1, s_data, 1);
3441
3442 target_data = data + (y * dest_line_length /*dest_stride * BPP/8*/);
3443
3444 for (int x = 0; x < target_width; x++) {
3445 unsigned int avgRed = 0;
3446 unsigned int avgGreen = 0;
3447 unsigned int avgBlue = 0;
3448 unsigned int pixel_count = 0;
3449 unsigned char *pix0 = s_data + BPP / 8 * ((int)(x * factor));
3450 y_offset = 0;
3451
3452 if ((x * Factor) < (Size_X - source.x)) {
3453 // determine average
3454 for (int y1 = 0; y1 < blur_factor; ++y1) {
3455 pixel = pix0 + (BPP / 8 * y_offset);
3456 for (int x1 = 0; x1 < blur_factor; ++x1) {
3457 avgRed += pixel[0];
3458 avgGreen += pixel[1];
3459 avgBlue += pixel[2];
3460
3461 pixel += BPP / 8;
3462
3463 pixel_count++;
3464 }
3465 y_offset += source.width;
3466 }
3467
3468 if (0 == pixel_count) // Protect
3469 pixel_count = 1;
3470
3471 target_data[0] = avgRed / pixel_count; // >> scounter;
3472 target_data[1] = avgGreen / pixel_count; // >> scounter;
3473 target_data[2] = avgBlue / pixel_count; // >> scounter;
3474 target_data += BPP / 8;
3475 } else {
3476 target_data[0] = 0;
3477 target_data[1] = 0;
3478 target_data[2] = 0;
3479 target_data += BPP / 8;
3480 }
3481
3482 } // for x
3483
3484 } // for y
3485
3486 } // SCALE_BILINEAR
3487
3488 else if (scale_type == RENDER_LODEF) {
3489 int get_bits_submap = 1;
3490
3491 int scaler = 16;
3492
3493 if (source.width > 32767) // High underscale can exceed signed math bits
3494 scaler = 8;
3495
3496 int wb_size = (Size_X) * ((/*Factor +*/ 1) * 2) * BPP / 8;
3497 s_data = (unsigned char *)malloc(wb_size); // work buffer
3498
3499 long x_delta = (source.width << scaler) / target_width;
3500 long y_delta = (source.height << scaler) / target_height;
3501
3502 int y = dest.y; // starting here
3503 long ys = dest.y * y_delta;
3504
3505 while (y < dest.y + dest.height) {
3506 // Read 1 line at the right place from the source
3507
3508 wxRect s1;
3509 s1.x = 0;
3510 s1.y = source.y + (ys >> scaler);
3511 s1.width = Size_X;
3512 s1.height = 1;
3513 GetChartBits(s1, s_data, get_bits_submap);
3514
3515 target_data = data + (y * dest_line_length /*dest_stride * BPP/8*/) +
3516 (dest.x * BPP / 8);
3517
3518 long x = (source.x << scaler) + (dest.x * x_delta);
3519 long sizex16 = Size_X << scaler;
3520 int xt = dest.x;
3521
3522 while ((xt < dest.x + dest.width) && (x < 0)) {
3523 target_data[0] = 0;
3524 target_data[1] = 0;
3525 target_data[2] = 0;
3526
3527 target_data += BPP / 8;
3528 x += x_delta;
3529 xt++;
3530 }
3531
3532 while ((xt < dest.x + dest.width) && (x < sizex16)) {
3533 unsigned char *src_pixel = &s_data[(x >> scaler) * BPP / 8];
3534
3535 target_data[0] = src_pixel[0];
3536 target_data[1] = src_pixel[1];
3537 target_data[2] = src_pixel[2];
3538
3539 target_data += BPP / 8;
3540 x += x_delta;
3541 xt++;
3542 }
3543
3544 while (xt < dest.x + dest.width) {
3545 target_data[0] = 0;
3546 target_data[1] = 0;
3547 target_data[2] = 0;
3548
3549 target_data += BPP / 8;
3550 xt++;
3551 }
3552
3553 y++;
3554 ys += y_delta;
3555 }
3556
3557 } // SCALE_SUBSAMP
3558
3559 } else // factor < 1, overzoom
3560 {
3561 int i = 0;
3562 int j = 0;
3563 unsigned char *target_line_start = NULL;
3564 unsigned char *target_data_x = NULL;
3565 int y_offset = 0;
3566
3567#ifdef __WXGTK__
3568 sigaction(SIGSEGV, NULL,
3569 &sa_all_previous); // save existing action for this signal
3570
3571 struct sigaction temp;
3572 sigaction(SIGSEGV, NULL, &temp); // inspect existing action for this signal
3573
3574 temp.sa_handler = catch_signals_chart; // point to my handler
3575 sigemptyset(&temp.sa_mask); // make the blocking set
3576 // empty, so that all
3577 // other signals will be
3578 // unblocked during my handler
3579 temp.sa_flags = 0;
3580 sigaction(SIGSEGV, &temp, NULL);
3581
3582 if (sigsetjmp(env_chart,
3583 1)) // Something in the below code block faulted....
3584 {
3585 sigaction(SIGSEGV, &sa_all_previous, NULL); // reset signal handler
3586
3587 wxString msg;
3588 msg.Printf(" Caught SIGSEGV on GetandScaleData, Factor < 1");
3589 wxLogMessage(msg);
3590
3591 msg.Printf(
3592 " m_raster_scale_factor: %g source.width: %d dest.y: %d "
3593 "dest.x: %d dest.width: %d dest.height: %d ",
3594 m_raster_scale_factor, source.width, dest.y, dest.x, dest.width,
3595 dest.height);
3596 wxLogMessage(msg);
3597
3598 msg.Printf(
3599 " i: %d j: %d dest_stride: %d target_line_start: %p "
3600 "target_data_x: %p y_offset: %d",
3601 i, j, dest_stride, target_line_start, target_data_x, y_offset);
3602 wxLogMessage(msg);
3603
3604 free(s_data);
3605 return true;
3606
3607 }
3608
3609 else
3610#endif
3611 {
3612
3613 double xd, yd;
3614 latlong_to_chartpix(m_vp_render_last.clat, m_vp_render_last.clon, xd, yd);
3615 double xrd =
3616 xd - (m_vp_render_last.pix_width * m_raster_scale_factor / 2);
3617 double yrd =
3618 yd - (m_vp_render_last.pix_height * m_raster_scale_factor / 2);
3619 double x_vernier = (xrd - wxRound(xrd));
3620 double y_vernier = (yrd - wxRound(yrd));
3621 int x_vernier_i = (int)wxRound(x_vernier / m_raster_scale_factor);
3622 int y_vernier_i = (int)wxRound(y_vernier / m_raster_scale_factor);
3623
3624 // Seems safe enough to read all the required data
3625 // Although we must adjust (increase) temporary allocation for negative
3626 // source.x and for vernier
3627 int sx = wxMax(source.x, 0);
3628 s_data = (unsigned char *)malloc((sx + source.width + 2) *
3629 (source.height + 2) * BPP / 8);
3630
3631 wxRect vsource = source;
3632 vsource.height += 2; // get more bits to allow for vernier
3633 vsource.width += 2;
3634 vsource.x -= 1;
3635 vsource.y -= 1;
3636
3637 GetChartBits(vsource, s_data, 1);
3638 unsigned char *source_data = s_data;
3639
3640 j = dest.y;
3641
3642 while (j < dest.y + dest.height) {
3643 y_offset = (int)((j - y_vernier_i) * m_raster_scale_factor) *
3644 vsource.width; // into the source data
3645
3646 target_line_start =
3647 target_data + (j * dest_line_length /*dest_stride * BPP / 8*/);
3648 target_data_x = target_line_start + ((dest.x) * BPP / 8);
3649
3650 i = dest.x;
3651
3652 // Check data bounds to be sure of not overrunning the upstream buffer
3653 if ((target_data_x + (dest.width * BPP / 8)) >
3654 (target_data + data_size)) {
3655 j = dest.y + dest.height;
3656 } else {
3657 while (i < dest.x + dest.width) {
3658 memcpy(target_data_x,
3659 source_data + BPP / 8 *
3660 (y_offset + (int)((i + x_vernier_i) *
3661 m_raster_scale_factor)),
3662 BPP / 8);
3663
3664 target_data_x += BPP / 8;
3665
3666 i++;
3667 }
3668 }
3669
3670 j++;
3671 }
3672 }
3673#ifdef __WXGTK__
3674 sigaction(SIGSEGV, &sa_all_previous, NULL); // reset signal handler
3675#endif
3676 }
3677
3678 free(s_data);
3679
3680 return true;
3681}
3682
3683bool ChartBaseBSB::GetChartBits(wxRect &source, unsigned char *pPix,
3684 int sub_samp) {
3685 wxCriticalSectionLocker locker(m_critSect);
3686
3687 int iy;
3688#define FILL_BYTE 0
3689
3690 // Decode the KAP file RLL stream into image pPix
3691
3692 unsigned char *pCP;
3693 pCP = pPix;
3694
3695 iy = source.y;
3696
3697 while (iy < source.y + source.height) {
3698 if ((iy >= 0) && (iy < Size_Y)) {
3699 if (source.x >= 0) {
3700 if ((source.x + source.width) > Size_X) {
3701 if ((Size_X - source.x) < 0)
3702 memset(pCP, FILL_BYTE, source.width * BPP / 8);
3703 else {
3704 BSBGetScanline(pCP, iy, source.x, Size_X, sub_samp);
3705 memset(pCP + (Size_X - source.x) * BPP / 8, FILL_BYTE,
3706 (source.x + source.width - Size_X) * BPP / 8);
3707 }
3708 } else
3709 BSBGetScanline(pCP, iy, source.x, source.x + source.width, sub_samp);
3710 } else {
3711 if ((source.width + source.x) >= 0) {
3712 // Special case, black on left side
3713 // must ensure that (black fill length % sub_samp) == 0
3714
3715 int xfill_corrected = -source.x + (source.x % sub_samp); //+ve
3716 memset(pCP, FILL_BYTE, (xfill_corrected * BPP / 8));
3717 BSBGetScanline(pCP + (xfill_corrected * BPP / 8), iy, 0,
3718 source.width + source.x, sub_samp);
3719
3720 } else {
3721 memset(pCP, FILL_BYTE, source.width * BPP / 8);
3722 }
3723 }
3724 }
3725
3726 else // requested y is off chart
3727 {
3728 memset(pCP, FILL_BYTE, source.width * BPP / 8);
3729 }
3730
3731 pCP += source.width * BPP / 8 * sub_samp;
3732
3733 iy += sub_samp;
3734 } // while iy
3735
3736 return true;
3737}
3738
3739//-----------------------------------------------------------------------------------------------
3740// BSB File Read Support
3741//-----------------------------------------------------------------------------------------------
3742
3743//-----------------------------------------------------------------------------------------------
3744// ReadBSBHdrLine
3745//
3746// Read and return count of a line of BSB header file
3747//-----------------------------------------------------------------------------------------------
3748
3749int ChartBaseBSB::ReadBSBHdrLine(wxInputStream *ifs, char *buf, int buf_len_max)
3750
3751{
3752 char read_char;
3753 char cr_test;
3754 int line_length = 0;
3755 char *lbuf;
3756
3757 lbuf = buf;
3758
3759 while (!ifs->Eof() && line_length < buf_len_max) {
3760 int c = ifs->GetC();
3761 if (c < 0) break;
3762 read_char = c;
3763 if (0x1A == read_char) {
3764 ifs->Ungetch(read_char);
3765 return (0);
3766 }
3767
3768 if (0 == read_char) // embedded erroneous unicode character?
3769 read_char = 0x20;
3770
3771 // Manage continued lines
3772 if (read_char == 10 || read_char == 13) {
3773 // Check to see if there is an extra CR
3774 cr_test = ifs->GetC();
3775 if (cr_test == 13) cr_test = ifs->GetC(); // skip any extra CR
3776
3777 if (cr_test != 10 && cr_test != 13) ifs->Ungetch(cr_test);
3778 read_char = '\n';
3779 }
3780
3781 // Look for continued lines, indicated by ' ' in first position
3782 if (read_char == '\n') {
3783 cr_test = 0;
3784 cr_test = ifs->GetC();
3785
3786 if (cr_test != ' ') {
3787 ifs->Ungetch(cr_test);
3788 *lbuf = '\0';
3789 return line_length;
3790 }
3791
3792 // Merge out leading spaces
3793 while (cr_test == ' ') cr_test = ifs->GetC();
3794 ifs->Ungetch(cr_test);
3795
3796 // Add a comma
3797 *lbuf = ',';
3798 lbuf++;
3799 }
3800
3801 else {
3802 *lbuf = read_char;
3803 lbuf++;
3804 line_length++;
3805 }
3806
3807 } // while
3808
3809 // Terminate line
3810 if (line_length) *(lbuf - 1) = '\0';
3811
3812 return line_length;
3813}
3814
3815//-----------------------------------------------------------------------
3816// Scan a BSB Scan Line from raw data
3817// Leaving stream pointer at start of next line
3818//-----------------------------------------------------------------------
3819int ChartBaseBSB::BSBScanScanline(wxInputStream *pinStream) {
3820 int nLineMarker, nValueShift, iPixel = 0;
3821 unsigned char byValueMask, byCountMask;
3822 unsigned char byNext;
3823 int coffset;
3824
3825 // if(1)
3826 {
3827 // Read the line number.
3828 nLineMarker = 0;
3829 do {
3830 byNext = pinStream->GetC();
3831 nLineMarker = nLineMarker * 128 + (byNext & 0x7f);
3832 } while ((byNext & 0x80) != 0);
3833
3834 // Setup masking values.
3835 nValueShift = 7 - nColorSize;
3836 byValueMask = (((1 << nColorSize)) - 1) << nValueShift;
3837 byCountMask = (1 << (7 - nColorSize)) - 1;
3838
3839 // Read and simulate expansion of runs.
3840
3841 while (((byNext = pinStream->GetC()) != 0) && (iPixel < Size_X)) {
3842 // int nPixValue;
3843 int nRunCount;
3844 // nPixValue = (byNext & byValueMask) >> nValueShift;
3845
3846 nRunCount = byNext & byCountMask;
3847
3848 while ((byNext & 0x80) != 0) {
3849 byNext = pinStream->GetC();
3850 nRunCount = nRunCount * 128 + (byNext & 0x7f);
3851 }
3852
3853 if (iPixel + nRunCount + 1 > Size_X) nRunCount = Size_X - iPixel - 1;
3854
3855 // Store nPixValue in the destination
3856 // memset(pCL, nPixValue, nRunCount+1);
3857 // pCL += nRunCount+1;
3858 iPixel += nRunCount + 1;
3859 }
3860 coffset = pinStream->TellI();
3861 }
3862
3863 return nLineMarker;
3864}
3865// MSVC compiler makes a bad decision about when to inline (or not) some
3866// intrinsics, like memset(). So,... Here is a little hand-crafted memset()
3867// substitue for known short strings. It will be inlined by MSVC compiler
3868// using /02 settings
3869
3870inline void memset_short(unsigned char *dst, unsigned char cbyte, int count) {
3871#if 0 // def __MSVC__
3872 __asm {
3873 pushf // save Direction flag
3874 cld // set direction "up"
3875 mov edi, dst
3876 mov ecx, count
3877 mov al, cbyte
3878 rep stosb
3879 popf
3880 }
3881#else
3882 memset(dst, cbyte, count);
3883#endif
3884}
3885// could use a larger value for slightly less ram but slower random access,
3886// this is chosen as it is also the opengl tile size so should work well
3887#define TILE_SIZE 512
3888
3889// #define USE_OLD_CACHE // removed this (and simplify code below) once the new
3890// method is verified #define PRINT_TIMINGS // enable for profiling
3891
3892#ifdef PRINT_TIMINGS
3893class OCPNStopWatch {
3894public:
3895 OCPNStopWatch() { Reset(); }
3896 void Reset() { clock_gettime(CLOCK_REALTIME, &tp); }
3897
3898 double Time() {
3899 timespec tp_end;
3900 clock_gettime(CLOCK_REALTIME, &tp_end);
3901 return (tp_end.tv_sec - tp.tv_sec) * 1.e3 +
3902 (tp_end.tv_nsec - tp.tv_nsec) / 1.e6;
3903 }
3904
3905private:
3906 timespec tp;
3907};
3908#endif
3909
3910#define FAIL \
3911 do { \
3912 free(pt->pTileOffset); \
3913 pt->pTileOffset = NULL; \
3914 free(pt->pPix); \
3915 pt->pPix = NULL; \
3916 pt->bValid = false; \
3917 return 0; \
3918 } while (0)
3919
3920//-----------------------------------------------------------------------
3921// Get a BSB Scan Line Using Cache and scan line index if available
3922//-----------------------------------------------------------------------
3923int ChartBaseBSB::BSBGetScanline(unsigned char *pLineBuf, int y, int xs, int xl,
3924 int sub_samp)
3925
3926{
3927 unsigned char *prgb = pLineBuf;
3928 int nValueShift;
3929 unsigned char byValueMask, byCountMask;
3930 unsigned char byNext;
3931 CachedLine *pt = NULL, cached_line;
3932 unsigned char *pCL;
3933 int rgbval;
3934 unsigned char *lp;
3935 int ix = xs;
3936 int pos = 0;
3937
3938 if (bUseLineCache && pLineCache) {
3939 // Is the requested line in the cache, and valid?
3940 pt = &pLineCache[y];
3941 } else {
3942 pt = &cached_line;
3943 pt->bValid = false;
3944 }
3945
3946#ifdef PRINT_TIMINGS
3947 OCPNStopWatch sw;
3948 static double ttime;
3949 static int cnt;
3950 cnt++;
3951#endif
3952
3953 if (!pt->bValid) // not valid, allocate
3954 {
3955 pt->size = pline_table[y + 1] - pline_table[y];
3956
3957#ifdef USE_OLD_CACHE
3958 pt->pPix = (unsigned char *)malloc(Size_X);
3959#else
3960 pt->pTileOffset = (TileOffsetCache *)calloc(
3961 sizeof(TileOffsetCache) * (Size_X / TILE_SIZE + 1), 1);
3962 pt->pPix = (unsigned char *)malloc(pt->size);
3963#endif
3964 if (pline_table[y] == 0 || pline_table[y + 1] == 0) FAIL;
3965
3966 // as of 2015, in wxWidgets buffered streams don't test for a zero seek
3967 // so we check here to possibly avoid this seek with a measured performance
3968 // gain
3969 if (ifs_bitmap->TellI() != pline_table[y] &&
3970 wxInvalidOffset == ifs_bitmap->SeekI(pline_table[y], wxFromStart))
3971 FAIL;
3972
3973#ifdef USE_OLD_CACHE
3974 if (pt->size > ifs_bufsize) {
3975 unsigned char *tmp = ifs_buf;
3976 if (!(ifs_buf = (unsigned char *)realloc(ifs_buf, pt->size))) {
3977 free(tmp);
3978 FAIL;
3979 }
3980 ifs_bufsize = pt->size;
3981 }
3982
3983 lp = ifs_buf;
3984#else
3985 lp = pt->pPix;
3986#endif
3987 ifs_bitmap->Read(lp, pt->size);
3988
3989#ifdef USE_OLD_CACHE
3990 pCL = pt->pPix;
3991#else
3992 if (!bUseLineCache) {
3993 ix = 0;
3994 // skip the line number.
3995 do {
3996 if (pos < pt->size) {
3997 byNext = *lp++;
3998 pos++;
3999 } else {
4000 byNext = 0;
4001 }
4002 } while ((byNext & 0x80) != 0);
4003 goto nocachestart;
4004 }
4005 pCL = ifs_buf;
4006
4007 if (Size_X > ifs_bufsize) {
4008 unsigned char *tmp = ifs_buf;
4009 if (!(ifs_buf = (unsigned char *)realloc(ifs_buf, Size_X))) {
4010 free(tmp);
4011 FAIL;
4012 }
4013 ifs_bufsize = Size_X;
4014 }
4015#endif
4016 // At this point, the unexpanded, raw line is at *lp, and the expansion
4017 // destination is pCL
4018
4019 // skip the line number.
4020 do {
4021 if (pos < pt->size) {
4022 byNext = *lp++;
4023 pos++;
4024 } else {
4025 byNext = 0;
4026 }
4027 } while ((byNext & 0x80) != 0);
4028
4029 // Setup masking values.
4030 nValueShift = 7 - nColorSize;
4031 byValueMask = (((1 << nColorSize)) - 1) << nValueShift;
4032 byCountMask = (1 << (7 - nColorSize)) - 1;
4033
4034 // Read and expand runs.
4035 unsigned int iPixel = 0;
4036
4037#ifndef USE_OLD_CACHE
4038 pt->pTileOffset[0].offset = lp - pt->pPix;
4039 pt->pTileOffset[0].pixel = 0;
4040 unsigned int tileindex = 1, nextTile = TILE_SIZE;
4041#endif
4042 unsigned int nRunCount;
4043 unsigned char *end = pt->pPix + pt->size;
4044 while (iPixel < (unsigned int)Size_X)
4045#ifdef USE_OLD_CACHE
4046 {
4047 nPixValue = (byNext & byValueMask) >> nValueShift;
4048
4049 nRunCount = byNext & byCountMask;
4050
4051 while ((byNext & 0x80) != 0) {
4052 byNext = *lp++;
4053 nRunCount = nRunCount * 128 + (byNext & 0x7f);
4054 }
4055
4056 nRunCount++;
4057
4058 if (iPixel + nRunCount >
4059 (unsigned int)Size_X) // protection against corrupt data
4060 nRunCount = nRunCount - iPixel;
4061
4062 // Store nPixValue in the destination
4063 memset_short(pCL + iPixel, nPixValue, nRunCount);
4064 iPixel += nRunCount;
4065 }
4066#else
4067 // build tile offset table for faster random access
4068 {
4069 if (pos < pt->size) {
4070 byNext = *lp++;
4071 pos++;
4072 } else {
4073 byNext = 0;
4074 }
4075 unsigned char *offset = lp - 1;
4076 if (byNext == 0 || lp == end) {
4077 // finished early...corrupt?
4078 while (tileindex < (unsigned int)Size_X / TILE_SIZE + 1) {
4079 pt->pTileOffset[tileindex].offset = pt->pTileOffset[0].offset;
4080 pt->pTileOffset[tileindex].pixel = 0;
4081 tileindex++;
4082 }
4083 break;
4084 }
4085
4086 nRunCount = byNext & byCountMask;
4087
4088 while ((byNext & 0x80) != 0) {
4089 if (pos < pt->size) {
4090 byNext = *lp++;
4091 pos++;
4092 } else {
4093 byNext = 0;
4094 }
4095 nRunCount = nRunCount * 128 + (byNext & 0x7f);
4096 }
4097
4098 nRunCount++;
4099
4100 if (iPixel + nRunCount >
4101 (unsigned int)Size_X) // protection against corrupt data
4102 nRunCount = Size_X - iPixel;
4103
4104 while (iPixel + nRunCount > nextTile) {
4105 pt->pTileOffset[tileindex].offset = offset - pt->pPix;
4106 pt->pTileOffset[tileindex].pixel = iPixel;
4107 tileindex++;
4108 nextTile += TILE_SIZE;
4109 }
4110 iPixel += nRunCount;
4111 }
4112#endif
4113
4114 pt->bValid = true;
4115 }
4116
4117 // Line is valid, de-reference thru proper pallete directly to target
4118
4119 if (xl > Size_X) xl = Size_X;
4120
4121#ifdef USE_OLD_CACHE
4122 pCL = pt->pPix + xs;
4123
4124 // Optimization for most usual case
4125 if ((BPP == 24) && (1 == sub_samp)) {
4126 ix = xs;
4127 while (ix < xl - 1) {
4128 unsigned char cur_by = *pCL;
4129 rgbval = (int)(pPalette[cur_by]);
4130 while ((ix < xl - 1)) {
4131 if (cur_by != *pCL) break;
4132 *((int *)prgb) = rgbval;
4133 prgb += 3;
4134 pCL++;
4135 ix++;
4136 }
4137 }
4138 } else {
4139 int dest_inc_val_bytes = (BPP / 8) * sub_samp;
4140 ix = xs;
4141 while (ix < xl - 1) {
4142 unsigned char cur_by = *pCL;
4143 rgbval = (int)(pPalette[cur_by]);
4144 while ((ix < xl - 1)) {
4145 if (cur_by != *pCL) break;
4146 *((int *)prgb) = rgbval;
4147 prgb += dest_inc_val_bytes;
4148 pCL += sub_samp;
4149 ix += sub_samp;
4150 }
4151 }
4152 }
4153
4154 // Get the last pixel explicitely
4155 // irrespective of the sub_sampling factor
4156
4157 if (xs < xl - 1) {
4158 unsigned char *pCLast = pt->pPix + (xl - 1);
4159 unsigned char *prgb_last = pLineBuf + ((xl - 1) - xs) * BPP / 8;
4160
4161 rgbval = (int)(pPalette[*pCLast]); // last pixel
4162 unsigned char a = rgbval & 0xff;
4163 *prgb_last++ = a;
4164 a = (rgbval >> 8) & 0xff;
4165 *prgb_last++ = a;
4166 a = (rgbval >> 16) & 0xff;
4167 *prgb_last = a;
4168 }
4169#else
4170 {
4171 int tileindex = xs / TILE_SIZE;
4172 int tileoffset = pt->pTileOffset[tileindex].offset;
4173
4174 lp = pt->pPix + tileoffset;
4175 pos = tileoffset;
4176 ix = pt->pTileOffset[tileindex].pixel;
4177 }
4178
4179nocachestart:
4180 nValueShift = 7 - nColorSize;
4181 byValueMask = (((1 << nColorSize)) - 1) << nValueShift;
4182 byCountMask = (1 << (7 - nColorSize)) - 1;
4183 int nPixValue = 0; // satisfy stupid compiler warning
4184 bool bLastPixValueValid = false;
4185 while (ix < xl - 1) {
4186 if (pos < pt->size) {
4187 byNext = *lp++;
4188 pos++;
4189 } else {
4190 break;
4191 }
4192
4193 nPixValue = (byNext & byValueMask) >> nValueShift;
4194 unsigned int nRunCount;
4195
4196 if (byNext == 0)
4197 nRunCount = xl - ix; // corrupted chart, just run to the end
4198 else {
4199 nRunCount = byNext & byCountMask;
4200 while ((byNext & 0x80) != 0) {
4201 if (pos < pt->size) {
4202 byNext = *lp++;
4203 pos++;
4204 } else {
4205 nRunCount = xl - ix; // corrupted chart, just run to the end
4206 break;
4207 }
4208 nRunCount = nRunCount * 128 + (byNext & 0x7f);
4209 }
4210
4211 nRunCount++;
4212 }
4213
4214 if (ix < xs) {
4215 if (ix + nRunCount <= (unsigned int)xs) {
4216 ix += nRunCount;
4217 continue;
4218 }
4219 nRunCount -= xs - ix;
4220 ix = xs;
4221 }
4222
4223 if (ix + nRunCount >= (unsigned int)xl) {
4224 nRunCount = xl - 1 - ix;
4225 bLastPixValueValid = true;
4226 }
4227
4228 rgbval = (int)(pPalette[nPixValue]);
4229
4230 // Optimization for most usual case
4231 // currently this is the only case possible...
4232 // if((BPP == 24) && (1 == sub_samp))
4233 {
4234 int count = nRunCount;
4235 if (count < 16) {
4236 // for short runs, use simple loop
4237 while (count--) {
4238 *(uint32_t *)prgb = rgbval;
4239 prgb += 3;
4240 }
4241 } else if (rgbval == 0 || rgbval == 0xffffff) {
4242 // optimization for black or white (could work for any gray too)
4243 memset(prgb, rgbval, nRunCount * 3);
4244 prgb += nRunCount * 3;
4245 } else {
4246 // note: this may not be optimal for all processors and compilers
4247 // I optimized for x86_64 using gcc with -O3
4248 // it is probably possible to gain even faster performance by ensuring
4249 // alignment to 16 or 32byte boundary (depending on processor) then
4250 // using inline assembly
4251
4252#ifdef __ARM_ARCH
4253 // ARM needs 8 byte alignment for *(uint64_T *x) = *(uint64_T *y)
4254 // because the compiler will (probably) use the ldrd/strd instuction
4255 // pair. So, advance the prgb pointer until it is 8-byte aligned, and
4256 // then carry on if enough bytes are left to process as 64 bit elements
4257
4258 if ((long)prgb & 7) {
4259 while (count--) {
4260 *(uint32_t *)prgb = rgbval;
4261 prgb += 3;
4262 if (!((long)prgb & 7)) {
4263 if (count >= 8) break;
4264 }
4265 }
4266 }
4267#endif
4268
4269 // fill first 24 bytes
4270 uint64_t *b = (uint64_t *)prgb;
4271 for (int i = 0; i < 8; i++) {
4272 *(uint32_t *)prgb = rgbval;
4273 prgb += 3;
4274 }
4275 count -= 8;
4276
4277 // fill in blocks of 24 bytes
4278 uint64_t *y = (uint64_t *)prgb;
4279 int count_d8 = count >> 3;
4280 prgb += 24 * count_d8;
4281 while (count_d8--) {
4282 *y++ = b[0];
4283 *y++ = b[1];
4284 *y++ = b[2];
4285 }
4286
4287 // fill remaining bytes
4288 int rcount = count & 0x7;
4289 while (rcount--) {
4290 *(uint32_t *)prgb = rgbval;
4291 prgb += 3;
4292 }
4293 }
4294 }
4295
4296 ix += nRunCount;
4297 }
4298
4299 // Get the last pixel explicitely
4300 // irrespective of the sub_sampling factor
4301
4302 if (ix < xl) {
4303 if (!bLastPixValueValid) {
4304 if (pos < pt->size) {
4305 byNext = *lp++;
4306 pos++;
4307 } else {
4308 byNext = 0;
4309 }
4310 nPixValue = (byNext & byValueMask) >> nValueShift;
4311 }
4312 rgbval = (int)(pPalette[nPixValue]); // last pixel
4313 unsigned char a = rgbval & 0xff;
4314
4315 *prgb++ = a;
4316 a = (rgbval >> 8) & 0xff;
4317 *prgb++ = a;
4318 a = (rgbval >> 16) & 0xff;
4319 *prgb = a;
4320 }
4321#endif
4322
4323#ifdef PRINT_TIMINGS
4324 ttime += sw.Time();
4325
4326 if (cnt == 500000) {
4327 static int d;
4328 printf("cache time: %d %f\n", d, ttime / 1000.0);
4329 cnt = 0;
4330 d++;
4331 // ttime = 0;
4332 }
4333#endif
4334
4335 if (!bUseLineCache) {
4336#ifndef USE_OLD_CACHE
4337 free(pt->pTileOffset);
4338#endif
4339 free(pt->pPix);
4340 }
4341
4342 return 1;
4343}
4344
4345int *ChartBaseBSB::GetPalettePtr(BSB_Color_Capability color_index) {
4346 if (pPalettes[color_index]) {
4347 if (palette_direction == PaletteFwd)
4348 return (int *)(pPalettes[color_index]->FwdPalette);
4349 else
4350 return (int *)(pPalettes[color_index]->RevPalette);
4351 } else
4352 return NULL;
4353}
4354
4355PaletteDir ChartBaseBSB::GetPaletteDir(void) {
4356 // make a pixel cache
4357 PixelCache *pc = new PixelCache(4, 4, BPP);
4358 RGBO r = pc->GetRGBO();
4359 delete pc;
4360
4361 if (r == RGB)
4362 return PaletteFwd;
4363 else
4364 return PaletteRev;
4365}
4366
4367bool ChartBaseBSB::AnalyzeSkew(void) {
4368 double lonmin = 1000;
4369 double lonmax = -1000;
4370 double latmin = 90.;
4371 double latmax = -90.;
4372
4373 int nlonmin, nlonmax;
4374 nlonmin = 0;
4375 nlonmax = 0;
4376
4377 if (0 == nRefpoint) // bad chart georef...
4378 return (1);
4379
4380 for (int n = 0; n < nRefpoint; n++) {
4381 // Longitude
4382 if (pRefTable[n].lonr > lonmax) {
4383 lonmax = pRefTable[n].lonr;
4384 nlonmax = n;
4385 }
4386 if (pRefTable[n].lonr < lonmin) {
4387 lonmin = pRefTable[n].lonr;
4388 nlonmin = n;
4389 }
4390
4391 // Latitude
4392 if (pRefTable[n].latr < latmin) {
4393 latmin = pRefTable[n].latr;
4394 }
4395 if (pRefTable[n].latr > latmax) {
4396 latmax = pRefTable[n].latr;
4397 }
4398 }
4399
4400 // Special case for charts which cross the IDL
4401 if ((lonmin * lonmax) < 0) {
4402 if (pRefTable[nlonmin].xr > pRefTable[nlonmax].xr) {
4403 // walk the reference table and add 360 to any longitude which is < 0
4404 for (int n = 0; n < nRefpoint; n++) {
4405 if (pRefTable[n].lonr < 0.0) pRefTable[n].lonr += 360.;
4406 }
4407
4408 // And recalculate the min/max
4409 lonmin = 1000;
4410 lonmax = -1000;
4411
4412 for (int n = 0; n < nRefpoint; n++) {
4413 // Longitude
4414 if (pRefTable[n].lonr > lonmax) {
4415 lonmax = pRefTable[n].lonr;
4416 nlonmax = n;
4417 }
4418 if (pRefTable[n].lonr < lonmin) {
4419 lonmin = pRefTable[n].lonr;
4420 nlonmin = n;
4421 }
4422
4423 // Latitude
4424 if (pRefTable[n].latr < latmin) {
4425 latmin = pRefTable[n].latr;
4426 }
4427 if (pRefTable[n].latr > latmax) {
4428 latmax = pRefTable[n].latr;
4429 }
4430 }
4431 }
4432 }
4433
4434 // Find the two REF points that are farthest apart
4435 double dist_max = 0.;
4436 int imax = 0;
4437 int jmax = 0;
4438
4439 for (int i = 0; i < nRefpoint; i++) {
4440 for (int j = i + 1; j < nRefpoint; j++) {
4441 double dx = pRefTable[i].xr - pRefTable[j].xr;
4442 double dy = pRefTable[i].yr - pRefTable[j].yr;
4443 double dist = (dx * dx) + (dy * dy);
4444 if (dist > dist_max) {
4445 dist_max = dist;
4446 imax = i;
4447 jmax = j;
4448 }
4449 }
4450 }
4451
4452 double apparent_skew = 0;
4453
4454 if (m_projection == PROJECTION_MERCATOR) {
4455 double easting0, easting1, northing0, northing1;
4456 // Get the Merc projection of the two REF points
4457 toSM_ECC(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4458 &easting0, &northing0);
4459 toSM_ECC(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4460 &easting1, &northing1);
4461
4462 double skew_proj =
4463 atan2((easting1 - easting0), (northing1 - northing0)) * 180. / PI;
4464 double skew_points = atan2((pRefTable[jmax].yr - pRefTable[imax].yr),
4465 (pRefTable[jmax].xr - pRefTable[imax].xr)) *
4466 180. / PI;
4467
4468 apparent_skew = skew_points - skew_proj + 90.;
4469
4470 // normalize to +/- 180.
4471 if (fabs(apparent_skew) > 180.) {
4472 if (apparent_skew < 0.)
4473 apparent_skew += 360.;
4474 else
4475 apparent_skew -= 360.;
4476 }
4477 }
4478
4479 else if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
4480 double easting0, easting1, northing0, northing1;
4481 // Get the TMerc projection of the two REF points
4482 toTM(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4483 &easting0, &northing0);
4484 toTM(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4485 &easting1, &northing1);
4486
4487 double skew_proj =
4488 atan2((easting1 - easting0), (northing1 - northing0)) * 180. / PI;
4489 double skew_points = atan2((pRefTable[jmax].yr - pRefTable[imax].yr),
4490 (pRefTable[jmax].xr - pRefTable[imax].xr)) *
4491 180. / PI;
4492
4493 apparent_skew = skew_points - skew_proj + 90.;
4494
4495 // normalize to +/- 180.
4496 if (fabs(apparent_skew) > 180.) {
4497 if (apparent_skew < 0.)
4498 apparent_skew += 360.;
4499 else
4500 apparent_skew -= 360.;
4501 }
4502
4503 if (fabs(apparent_skew - m_Chart_Skew) > 2) { // measured skew OK?
4504 // it may be that the projection longitude is simply wrong.
4505 // Check it.
4506 double dskew = fabs(apparent_skew - m_Chart_Skew);
4507 if ((m_proj_lon < lonmin) || (m_proj_lon > lonmax)) {
4508 // Try a projection longitude that is mid-meridian in the chart.
4509 double tentative_proj_lon = (lonmin + lonmax) / 2.;
4510
4511 toTM(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat,
4512 tentative_proj_lon, &easting0, &northing0);
4513 toTM(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat,
4514 tentative_proj_lon, &easting1, &northing1);
4515
4516 skew_proj =
4517 atan2((easting1 - easting0), (northing1 - northing0)) * 180. / PI;
4518 skew_points = atan2((pRefTable[jmax].yr - pRefTable[imax].yr),
4519 (pRefTable[jmax].xr - pRefTable[imax].xr)) *
4520 180. / PI;
4521
4522 apparent_skew = skew_points - skew_proj + 90.;
4523
4524 // normalize to +/- 180.
4525 if (fabs(apparent_skew) > 180.) {
4526 if (apparent_skew < 0.)
4527 apparent_skew += 360.;
4528 else
4529 apparent_skew -= 360.;
4530 }
4531
4532 // Better? If so, adopt the adjusted projection longitude
4533 if (fabs(apparent_skew - m_Chart_Skew) < dskew) {
4534 m_proj_lon = tentative_proj_lon;
4535 }
4536 }
4537 }
4538 } else // For all other projections, assume that skew specified in header is
4539 // correct
4540 apparent_skew = m_Chart_Skew;
4541
4542 if (fabs(apparent_skew - m_Chart_Skew) >
4543 2) { // measured skew is more than 2 degrees
4544 m_Chart_Skew = apparent_skew; // different from stated skew
4545
4546 wxString msg = " Warning: Skew override on chart ";
4547 msg.Append(m_FullPath);
4548 wxString msg1;
4549 msg1.Printf(" is %5g degrees", apparent_skew);
4550 msg.Append(msg1);
4551
4552 wxLogMessage(msg);
4553
4554 return false;
4555 }
4556
4557 return true;
4558}
4559
4560int ChartBaseBSB::AnalyzeRefpoints(bool b_testSolution) {
4561 int i, n;
4562 double elt, elg;
4563
4564 // Calculate the max/min reference points
4565
4566 float lonmin = 1000;
4567 float lonmax = -1000;
4568 float latmin = 90.;
4569 float latmax = -90.;
4570
4571 int plonmin = 100000;
4572 int plonmax = 0;
4573 int platmin = 100000;
4574 int platmax = 0;
4575 int nlonmin = 0, nlonmax = 0;
4576
4577 if (0 == nRefpoint) // bad chart georef...
4578 return (1);
4579
4580 for (n = 0; n < nRefpoint; n++) {
4581 // Longitude
4582 if (pRefTable[n].lonr > lonmax) {
4583 lonmax = pRefTable[n].lonr;
4584 plonmax = (int)pRefTable[n].xr;
4585 nlonmax = n;
4586 }
4587 if (pRefTable[n].lonr < lonmin) {
4588 lonmin = pRefTable[n].lonr;
4589 plonmin = (int)pRefTable[n].xr;
4590 nlonmin = n;
4591 }
4592
4593 // Latitude
4594 if (pRefTable[n].latr < latmin) {
4595 latmin = pRefTable[n].latr;
4596 platmin = (int)pRefTable[n].yr;
4597 }
4598 if (pRefTable[n].latr > latmax) {
4599 latmax = pRefTable[n].latr;
4600 platmax = (int)pRefTable[n].yr;
4601 }
4602 }
4603
4604 // Special case for charts which cross the IDL
4605 if ((lonmin * lonmax) < 0) {
4606 if (pRefTable[nlonmin].xr > pRefTable[nlonmax].xr) {
4607 // walk the reference table and add 360 to any longitude which is < 0
4608 for (n = 0; n < nRefpoint; n++) {
4609 if (pRefTable[n].lonr < 0.0) pRefTable[n].lonr += 360.;
4610 }
4611
4612 // And recalculate the min/max
4613 lonmin = 1000;
4614 lonmax = -1000;
4615
4616 for (n = 0; n < nRefpoint; n++) {
4617 // Longitude
4618 if (pRefTable[n].lonr > lonmax) {
4619 lonmax = pRefTable[n].lonr;
4620 plonmax = (int)pRefTable[n].xr;
4621 nlonmax = n;
4622 }
4623 if (pRefTable[n].lonr < lonmin) {
4624 lonmin = pRefTable[n].lonr;
4625 plonmin = (int)pRefTable[n].xr;
4626 nlonmin = n;
4627 }
4628
4629 // Latitude
4630 if (pRefTable[n].latr < latmin) {
4631 latmin = pRefTable[n].latr;
4632 platmin = (int)pRefTable[n].yr;
4633 }
4634 if (pRefTable[n].latr > latmax) {
4635 latmax = pRefTable[n].latr;
4636 platmax = (int)pRefTable[n].yr;
4637 }
4638 }
4639 }
4640 }
4641
4642 // Build the Control Point Structure, etc
4643 cPoints.count = nRefpoint;
4644 if (cPoints.status) {
4645 // AnalyzeRefpoints can be called twice
4646 free(cPoints.tx);
4647 free(cPoints.ty);
4648 free(cPoints.lon);
4649 free(cPoints.lat);
4650
4651 free(cPoints.pwx);
4652 free(cPoints.wpx);
4653 free(cPoints.pwy);
4654 free(cPoints.wpy);
4655 }
4656
4657 cPoints.tx = (double *)malloc(nRefpoint * sizeof(double));
4658 cPoints.ty = (double *)malloc(nRefpoint * sizeof(double));
4659 cPoints.lon = (double *)malloc(nRefpoint * sizeof(double));
4660 cPoints.lat = (double *)malloc(nRefpoint * sizeof(double));
4661
4662 cPoints.pwx = (double *)malloc(12 * sizeof(double));
4663 cPoints.wpx = (double *)malloc(12 * sizeof(double));
4664 cPoints.pwy = (double *)malloc(12 * sizeof(double));
4665 cPoints.wpy = (double *)malloc(12 * sizeof(double));
4666 cPoints.status = 1;
4667
4668 // Find the two REF points that are farthest apart
4669 double dist_max = 0.;
4670 int imax = 0;
4671 int jmax = 0;
4672
4673 for (i = 0; i < nRefpoint; i++) {
4674 for (int j = i + 1; j < nRefpoint; j++) {
4675 double dx = pRefTable[i].xr - pRefTable[j].xr;
4676 double dy = pRefTable[i].yr - pRefTable[j].yr;
4677 double dist = (dx * dx) + (dy * dy);
4678 if (dist > dist_max) {
4679 dist_max = dist;
4680 imax = i;
4681 jmax = j;
4682 }
4683 }
4684 }
4685
4686 // Georef solution depends on projection type
4687
4688 if (m_projection == PROJECTION_TRANSVERSE_MERCATOR) {
4689 double easting0, easting1, northing0, northing1;
4690 // Get the TMerc projection of the two REF points
4691 toTM(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4692 &easting0, &northing0);
4693 toTM(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4694 &easting1, &northing1);
4695
4696 // Calculate the scale factor using exact REF point math
4697 double dx2 = (pRefTable[jmax].xr - pRefTable[imax].xr) *
4698 (pRefTable[jmax].xr - pRefTable[imax].xr);
4699 double dy2 = (pRefTable[jmax].yr - pRefTable[imax].yr) *
4700 (pRefTable[jmax].yr - pRefTable[imax].yr);
4701 double dn2 = (northing1 - northing0) * (northing1 - northing0);
4702 double de2 = (easting1 - easting0) * (easting1 - easting0);
4703
4704 m_ppm_avg = sqrt(dx2 + dy2) / sqrt(dn2 + de2);
4705
4706 // Set up and solve polynomial solution for pix<->east/north as projected
4707 // Fill the cpoints structure with pixel points and transformed lat/lon
4708
4709 for (int n = 0; n < nRefpoint; n++) {
4710 double easting, northing;
4711 toTM(pRefTable[n].latr, pRefTable[n].lonr, m_proj_lat, m_proj_lon,
4712 &easting, &northing);
4713
4714 cPoints.tx[n] = pRefTable[n].xr;
4715 cPoints.ty[n] = pRefTable[n].yr;
4716 cPoints.lon[n] = easting;
4717 cPoints.lat[n] = northing;
4718 }
4719
4720 // Helper parameters
4721 cPoints.txmax = plonmax;
4722 cPoints.txmin = plonmin;
4723 cPoints.tymax = platmax;
4724 cPoints.tymin = platmin;
4725 toTM(latmax, lonmax, m_proj_lat, m_proj_lon, &cPoints.lonmax,
4726 &cPoints.latmax);
4727 toTM(latmin, lonmin, m_proj_lat, m_proj_lon, &cPoints.lonmin,
4728 &cPoints.latmin);
4729
4730 Georef_Calculate_Coefficients_Proj(&cPoints);
4731
4732 }
4733
4734 else if (m_projection == PROJECTION_MERCATOR) {
4735 double easting0, easting1, northing0, northing1;
4736 // Get the Merc projection of the two REF points
4737 toSM_ECC(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4738 &easting0, &northing0);
4739 toSM_ECC(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4740 &easting1, &northing1);
4741
4742 // Calculate the scale factor using exact REF point math
4743 // double dx = (pRefTable[jmax].xr - pRefTable[imax].xr);
4744 // double de = (easting1 - easting0);
4745 // m_ppm_avg = fabs(dx / de);
4746
4747 double dx2 = (pRefTable[jmax].xr - pRefTable[imax].xr) *
4748 (pRefTable[jmax].xr - pRefTable[imax].xr);
4749 double dy2 = (pRefTable[jmax].yr - pRefTable[imax].yr) *
4750 (pRefTable[jmax].yr - pRefTable[imax].yr);
4751 double dn2 = (northing1 - northing0) * (northing1 - northing0);
4752 double de2 = (easting1 - easting0) * (easting1 - easting0);
4753
4754 m_ppm_avg = sqrt(dx2 + dy2) / sqrt(dn2 + de2);
4755
4756 // Set up and solve polynomial solution for pix<->east/north as projected
4757 // Fill the cpoints structure with pixel points and transformed lat/lon
4758
4759 for (int n = 0; n < nRefpoint; n++) {
4760 double easting, northing;
4761 toSM_ECC(pRefTable[n].latr, pRefTable[n].lonr, m_proj_lat, m_proj_lon,
4762 &easting, &northing);
4763
4764 cPoints.tx[n] = pRefTable[n].xr;
4765 cPoints.ty[n] = pRefTable[n].yr;
4766 cPoints.lon[n] = easting;
4767 cPoints.lat[n] = northing;
4768 // printf(" x: %g y: %g east: %g north:
4769 // %g\n",pRefTable[n].xr, pRefTable[n].yr, easting,
4770 // northing);
4771 }
4772
4773 // Helper parameters
4774 cPoints.txmax = plonmax;
4775 cPoints.txmin = plonmin;
4776 cPoints.tymax = platmax;
4777 cPoints.tymin = platmin;
4778 toSM_ECC(latmax, lonmax, m_proj_lat, m_proj_lon, &cPoints.lonmax,
4779 &cPoints.latmax);
4780 toSM_ECC(latmin, lonmin, m_proj_lat, m_proj_lon, &cPoints.lonmin,
4781 &cPoints.latmin);
4782
4783 Georef_Calculate_Coefficients_Proj(&cPoints);
4784
4785 // for(int h=0 ; h < 10 ; h++)
4786 // printf("pix to east %d %g\n", h, cPoints.pwx[h]); //
4787 // pix to lon
4788 // for(int h=0 ; h < 10 ; h++)
4789 // printf("east to pix %d %g\n", h, cPoints.wpx[h]); //
4790 // lon to pix
4791
4792 /*
4793 if ((0 != m_Chart_DU ) && (0 != m_Chart_Scale))
4794 {
4795 double m_ppm_avg1 = m_Chart_DU * 39.37 / m_Chart_Scale;
4796 m_ppm_avg1 *= cos(m_proj_lat * PI / 180.); // correct to
4797 chart centroid
4798
4799 printf("BSB chart ppm_avg:%g ppm_avg1:%g\n", m_ppm_avg,
4800 m_ppm_avg1); m_ppm_avg = m_ppm_avg1;
4801 }
4802 */
4803 }
4804
4805 else if (m_projection == PROJECTION_POLYCONIC) {
4806 // This is interesting
4807 // On some BSB V 1.0 Polyconic charts (e.g. 14500_1, 1995), the projection
4808 // parameter Which is taken to be the central meridian of the projection
4809 // is of the wrong sign....
4810
4811 // We check for this case, and make a correction if necessary.....
4812 // Obviously, the projection meridian should be on the chart, i.e. between
4813 // the min and max longitudes....
4814 double proj_meridian = m_proj_lon;
4815
4816 if ((pRefTable[nlonmax].lonr >= -proj_meridian) &&
4817 (-proj_meridian >= pRefTable[nlonmin].lonr))
4818 m_proj_lon = -m_proj_lon;
4819
4820 double easting0, easting1, northing0, northing1;
4821 // Get the Poly projection of the two REF points
4822 toPOLY(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4823 &easting0, &northing0);
4824 toPOLY(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4825 &easting1, &northing1);
4826
4827 // Calculate the scale factor using exact REF point math
4828 double dx2 = (pRefTable[jmax].xr - pRefTable[imax].xr) *
4829 (pRefTable[jmax].xr - pRefTable[imax].xr);
4830 double dy2 = (pRefTable[jmax].yr - pRefTable[imax].yr) *
4831 (pRefTable[jmax].yr - pRefTable[imax].yr);
4832 double dn2 = (northing1 - northing0) * (northing1 - northing0);
4833 double de2 = (easting1 - easting0) * (easting1 - easting0);
4834
4835 m_ppm_avg = sqrt(dx2 + dy2) / sqrt(dn2 + de2);
4836
4837 // Sanity check
4838 // double ref_dist = DistGreatCircle(pRefTable[imax].latr,
4839 // pRefTable[imax].lonr, pRefTable[jmax].latr,
4840 // pRefTable[jmax].lonr); ref_dist *= 1852; //To Meters double
4841 // ref_dist_transform = sqrt(dn2 + de2); //Also meters
4842 // double error = (ref_dist - ref_dist_transform)/ref_dist;
4843
4844 // Set up and solve polynomial solution for pix<->cartesian east/north as
4845 // projected
4846 // Fill the cpoints structure with pixel points and transformed lat/lon
4847
4848 for (int n = 0; n < nRefpoint; n++) {
4849 double easting, northing;
4850 toPOLY(pRefTable[n].latr, pRefTable[n].lonr, m_proj_lat, m_proj_lon,
4851 &easting, &northing);
4852
4853 // Round trip check for debugging....
4854 // double lat, lon;
4855 // fromPOLY(easting, northing, m_proj_lat, m_proj_lon,
4856 // &lat, &lon);
4857
4858 cPoints.tx[n] = pRefTable[n].xr;
4859 cPoints.ty[n] = pRefTable[n].yr;
4860 cPoints.lon[n] = easting;
4861 cPoints.lat[n] = northing;
4862 // printf(" x: %g y: %g east: %g north:
4863 // %g\n",pRefTable[n].xr, pRefTable[n].yr, easting,
4864 // northing);
4865 }
4866
4867 // Helper parameters
4868 cPoints.txmax = plonmax;
4869 cPoints.txmin = plonmin;
4870 cPoints.tymax = platmax;
4871 cPoints.tymin = platmin;
4872 toPOLY(latmax, lonmax, m_proj_lat, m_proj_lon, &cPoints.lonmax,
4873 &cPoints.latmax);
4874 toPOLY(latmin, lonmin, m_proj_lat, m_proj_lon, &cPoints.lonmin,
4875 &cPoints.latmin);
4876
4877 Georef_Calculate_Coefficients_Proj(&cPoints);
4878
4879 // for(int h=0 ; h < 10 ; h++)
4880 // printf("pix to east %d %g\n", h, cPoints.pwx[h]); //
4881 // pix to lon
4882 // for(int h=0 ; h < 10 ; h++)
4883 // printf("east to pix %d %g\n", h, cPoints.wpx[h]); //
4884 // lon to pix
4885
4886 }
4887
4888 // Any other projection had better have embedded coefficients
4889 else if (bHaveEmbeddedGeoref) {
4890 // Use a Mercator Projection to get a rough ppm.
4891 double easting0, easting1, northing0, northing1;
4892 // Get the Merc projection of the two REF points
4893 toSM_ECC(pRefTable[imax].latr, pRefTable[imax].lonr, m_proj_lat, m_proj_lon,
4894 &easting0, &northing0);
4895 toSM_ECC(pRefTable[jmax].latr, pRefTable[jmax].lonr, m_proj_lat, m_proj_lon,
4896 &easting1, &northing1);
4897
4898 // Calculate the scale factor using exact REF point math in x(longitude)
4899 // direction
4900
4901 double dx = (pRefTable[jmax].xr - pRefTable[imax].xr);
4902 double de = (easting1 - easting0);
4903
4904 m_ppm_avg = fabs(dx / de);
4905
4906 m_ExtraInfo = "---<<< Warning: Chart georef accuracy may be poor. >>>---";
4907 }
4908
4909 else
4910 m_ppm_avg = 1.0; // absolute fallback to prevent div-0 errors
4911
4912#if 0
4913 // Alternate Skew verification
4914 ViewPort vps;
4915 vps.clat = pRefTable[0].latr;
4916 vps.clon = pRefTable[0].lonr;
4917 vps.view_scale_ppm = m_ppm_avg;
4918 vps.skew = 0.;
4919 vps.pix_width = 1000;
4920 vps.pix_height = 1000;
4921
4922 int x1, y1, x2, y2;
4923 latlong_to_pix_vp(latmin, (lonmax + lonmin)/2., x1, y1, vps);
4924 latlong_to_pix_vp(latmax, (lonmax + lonmin)/2., x2, y2, vps);
4925
4926 double apparent_skew = (atan2( (y2-y1), (x2-x1) ) * 180./PI) + 90.;
4927 if(apparent_skew < 0.)
4928 apparent_skew += 360;
4929 if(apparent_skew > 360.)
4930 apparent_skew -= 360;
4931
4932 if(fabs( apparent_skew - m_Chart_Skew ) > 2) { // measured skew is more than 2 degress different
4933 m_Chart_Skew = apparent_skew;
4934 }
4935#endif
4936
4937 if (!b_testSolution) return (0);
4938
4939 // Do a last little test using a synthetic ViewPort of nominal size.....
4940 ViewPort vp;
4941 vp.clat = pRefTable[0].latr;
4942 vp.clon = pRefTable[0].lonr;
4943 vp.view_scale_ppm = m_ppm_avg;
4944 vp.skew = 0.;
4945 vp.pix_width = 1000;
4946 vp.pix_height = 1000;
4947 // vp.rv_rect = wxRect(0,0, vp.pix_width, vp.pix_height);
4948 SetVPRasterParms(vp);
4949
4950 double xpl_err_max = 0;
4951 double ypl_err_max = 0;
4952 int px, py;
4953
4954 int pxref, pyref;
4955 pxref = (int)pRefTable[0].xr;
4956 pyref = (int)pRefTable[0].yr;
4957
4958 for (i = 0; i < nRefpoint; i++) {
4959 px = (int)(vp.pix_width / 2 + pRefTable[i].xr) - pxref;
4960 py = (int)(vp.pix_height / 2 + pRefTable[i].yr) - pyref;
4961
4962 vp_pix_to_latlong(vp, px, py, &elt, &elg);
4963
4964 double lat_error = elt - pRefTable[i].latr;
4965 pRefTable[i].ypl_error = lat_error;
4966
4967 double lon_error = elg - pRefTable[i].lonr;
4968
4969 // Longitude errors could be compounded by prior adjustment to ref points
4970 if (fabs(lon_error) > 180.) {
4971 if (lon_error > 0.)
4972 lon_error -= 360.;
4973 else if (lon_error < 0.)
4974 lon_error += 360.;
4975 }
4976 pRefTable[i].xpl_error = lon_error;
4977
4978 if (fabs(pRefTable[i].ypl_error) > fabs(ypl_err_max))
4979 ypl_err_max = pRefTable[i].ypl_error;
4980 if (fabs(pRefTable[i].xpl_error) > fabs(xpl_err_max))
4981 xpl_err_max = pRefTable[i].xpl_error;
4982 }
4983
4984 Chart_Error_Factor = fmax(fabs(xpl_err_max / (lonmax - lonmin)),
4985 fabs(ypl_err_max / (latmax - latmin)));
4986 double chart_error_meters =
4987 fmax(fabs(xpl_err_max * 60. * 1852.), fabs(ypl_err_max * 60. * 1852.));
4988 // calculate a nominal pixel error
4989 // Assume a modern display has about 4000 pixels/meter.
4990 // Assume the chart is to be displayed at nominal printed scale
4991 double chart_error_pixels = chart_error_meters * 4000. / m_Chart_Scale;
4992
4993 // Good enough for navigation?
4994 int max_pixel_error = 4;
4995
4996 if (chart_error_pixels > max_pixel_error) {
4997 wxString msg =
4998 " VP Final Check: Georeference Chart_Error_Factor on chart ";
4999 msg.Append(m_FullPath);
5000 wxString msg1;
5001 msg1.Printf(" is %5g \n nominal pixel error is: %5g",
5002 Chart_Error_Factor, chart_error_pixels);
5003 msg.Append(msg1);
5004
5005 wxLogMessage(msg);
5006
5007 m_ExtraInfo = "---<<< Warning: Chart georef accuracy is poor. >>>---";
5008 }
5009
5010 // Try again with my calculated georef
5011 // This problem was found on NOAA 514_1.KAP. The embedded coefficients are
5012 // just wrong....
5013 if ((chart_error_pixels > max_pixel_error) && bHaveEmbeddedGeoref) {
5014 wxString msg =
5015 " Trying again with internally calculated georef solution ";
5016 wxLogMessage(msg);
5017
5018 bHaveEmbeddedGeoref = false;
5019 SetVPRasterParms(vp);
5020
5021 xpl_err_max = 0;
5022 ypl_err_max = 0;
5023
5024 pxref = (int)pRefTable[0].xr;
5025 pyref = (int)pRefTable[0].yr;
5026
5027 for (i = 0; i < nRefpoint; i++) {
5028 px = (int)(vp.pix_width / 2 + pRefTable[i].xr) - pxref;
5029 py = (int)(vp.pix_height / 2 + pRefTable[i].yr) - pyref;
5030
5031 vp_pix_to_latlong(vp, px, py, &elt, &elg);
5032
5033 double lat_error = elt - pRefTable[i].latr;
5034 pRefTable[i].ypl_error = lat_error;
5035
5036 double lon_error = elg - pRefTable[i].lonr;
5037
5038 // Longitude errors could be compounded by prior adjustment to ref points
5039 if (fabs(lon_error) > 180.) {
5040 if (lon_error > 0.)
5041 lon_error -= 360.;
5042 else if (lon_error < 0.)
5043 lon_error += 360.;
5044 }
5045 pRefTable[i].xpl_error = lon_error;
5046
5047 if (fabs(pRefTable[i].ypl_error) > fabs(ypl_err_max))
5048 ypl_err_max = pRefTable[i].ypl_error;
5049 if (fabs(pRefTable[i].xpl_error) > fabs(xpl_err_max))
5050 xpl_err_max = pRefTable[i].xpl_error;
5051 }
5052
5053 Chart_Error_Factor = fmax(fabs(xpl_err_max / (lonmax - lonmin)),
5054 fabs(ypl_err_max / (latmax - latmin)));
5055
5056 chart_error_meters =
5057 fmax(fabs(xpl_err_max * 60. * 1852.), fabs(ypl_err_max * 60. * 1852.));
5058 chart_error_pixels = chart_error_meters * 4000. / m_Chart_Scale;
5059
5060 // Good enough for navigation?
5061 if (chart_error_pixels > max_pixel_error) {
5062 wxString msg =
5063 " VP Final Check with internal georef: Georeference "
5064 "Chart_Error_Factor on chart ";
5065 msg.Append(m_FullPath);
5066 wxString msg1;
5067 msg1.Printf(" is %5g\n nominal pixel error is: %5g",
5068 Chart_Error_Factor, chart_error_pixels);
5069 msg.Append(msg1);
5070
5071 wxLogMessage(msg);
5072
5073 m_ExtraInfo = "---<<< Warning: Chart georef accuracy is poor. >>>---";
5074 } else {
5075 wxString msg = " Result: OK, Internal georef solution used.";
5076
5077 wxLogMessage(msg);
5078
5079 m_ExtraInfo = "";
5080 }
5081 }
5082
5083 return (0);
5084}
5085
5086double ChartBaseBSB::AdjustLongitude(double lon) {
5087 double lond = (m_LonMin + m_LonMax) / 2 - lon;
5088 if (lond > 180)
5089 return lon + 360;
5090 else if (lond < -180)
5091 return lon - 360;
5092 return lon;
5093}
5094
5095/*
5096 * Extracted from bsb_io.c - implementation of libbsb reading and writing
5097 *
5098 * Copyright (C) 2000 Stuart Cunningham <stuart_hc@users.sourceforge.net>
5099 *
5100 * This library is free software; you can redistribute it and/or
5101 * modify it under the terms of the GNU Lesser General Public
5102 * License as published by the Free Software Foundation; either
5103 * version 2.1 of the License, or (at your option) any later version.
5104 *
5105 * This library is distributed in the hope that it will be useful,
5106 * but WITHOUT ANY WARRANTY; without even the implied warranty of
5107 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
5108 * Lesser General Public License for more details.
5109 *
5110 * You should have received a copy of the GNU Lesser General Public
5111 * License along with this library; if not, write to the Free Software
5112 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
5113 *
5114 *
5115 */
5116
5126static double polytrans(double *coeff, double lon, double lat) {
5127 double ret = coeff[0] + coeff[1] * lon + coeff[2] * lat;
5128 ret += coeff[3] * lon * lon;
5129 ret += coeff[4] * lon * lat;
5130 ret += coeff[5] * lat * lat;
5131 ret += coeff[6] * lon * lon * lon;
5132 ret += coeff[7] * lon * lon * lat;
5133 ret += coeff[8] * lon * lat * lat;
5134 ret += coeff[9] * lat * lat * lat;
5135 // ret += coeff[10]*lat*lat*lat*lat;
5136 // ret += coeff[11]*lat*lat*lat*lat*lat;
5137
5138 // for(int n = 0 ; n < 10 ; n++)
5139 // printf(" %g\n", coeff[n]);
5140
5141 return ret;
5142}
5143
5144#if 0
5156extern int bsb_LLtoXY(BSBImage *p, double lon, double lat, int* x, int* y)
5157{
5158 /* change longitude phase (CPH) */
5159 lon = (lon < 0) ? lon + p->cph : lon - p->cph;
5160 double xd = polytrans( p->wpx, lon, lat );
5161 double yd = polytrans( p->wpy, lon, lat );
5162 *x = (int)(xd + 0.5);
5163 *y = (int)(yd + 0.5);
5164 return 1;
5165}
5166
5178extern int bsb_XYtoLL(BSBImage *p, int x, int y, double* lonout, double* latout)
5179{
5180 double lon = polytrans( p->pwx, x, y );
5181 lon = (lon < 0) ? lon + p->cph : lon - p->cph;
5182 *lonout = lon;
5183 *latout = polytrans( p->pwy, x, y );
5184 return 1;
5185}
5186
5187#endif
XZ compressed charts support.
BSB chart management.
virtual double GetNearestPreferredScalePPM(double target_scale_ppm)
Find the nearest preferred viewport scale (in pixels/meter) for this chart.
double GetClosestValidNaturalScalePPM(double target_scale, double scale_factor_min, double scale_factor_max)
Find closest valid scale that's a power of 2 multiple of chart's native scale.
An iterator class for OCPNRegion.
Definition OCPNRegion.h:156
A wrapper class for wxRegion with additional functionality.
Definition OCPNRegion.h:56
ViewPort - Core geographic projection and coordinate transformation engine.
Definition viewport.h:84
double view_scale_ppm
Requested view scale in physical pixels per meter (ppm), before applying projections.
Definition viewport.h:232
int pix_height
Height of the viewport in physical pixels.
Definition viewport.h:261
int pix_width
Width of the viewport in physical pixels.
Definition viewport.h:259
double skew
Angular distortion (shear transform) applied to the viewport in radians.
Definition viewport.h:240
double clon
Center longitude of the viewport in degrees.
Definition viewport.h:227
double clat
Center latitude of the viewport in degrees.
Definition viewport.h:225
double lonmin
Minimum longitude in reference data.
Definition georef.h:89
int count
Number of reference points used.
Definition georef.h:72
double lonmax
Maximum longitude in reference data.
Definition georef.h:88
int txmin
Minimum x value in target space.
Definition georef.h:86
double * pwx
Polynomial coefficients for pixel-to-world longitude transformation.
Definition georef.h:78
double * tx
Array of x-coordinates in target (typically pixel) space.
Definition georef.h:74
double * wpx
Polynomial coefficients for world-to-pixel x transformation.
Definition georef.h:82
double * lat
Array of latitudes corresponding to reference points.
Definition georef.h:77
double * lon
Array of longitudes corresponding to reference points.
Definition georef.h:76
int status
Status of the georeferencing (0 = valid, other values indicate errors)
Definition georef.h:70
int tymin
Minimum y value in target space.
Definition georef.h:87
double * ty
Array of y-coordinates in target (typically pixel) space.
Definition georef.h:75
int tymax
Maximum y value in target space.
Definition georef.h:85
double latmin
Minimum latitude in reference data.
Definition georef.h:91
int txmax
Maximum x value in target space.
Definition georef.h:84
double latmax
Maximum latitude in reference data.
Definition georef.h:90
double * wpy
Polynomial coefficients for world-to-pixel y transformation.
Definition georef.h:83
double * pwy
Polynomial coefficients for pixel-to-world latitude transformation.
Definition georef.h:80