← Back to branch summary

~ubuntu-branches/debian/sid/gdal/sid

« back to all changes in this revision

Viewing changes to frmts/pcidsk/sdk/segment/cpcidskgcp2segment.cpp

  • Committer: Package Import Robot
  • Author(s): Francesco Paolo Lovergine
  • Date: 2012-05-07 15:04:42 UTC
  • mfrom: (5.5.16 experimental)
  • Revision ID: package-import@ubuntu.com-20120507150442-2eks97loeh6rq005
Tags: 1.9.0-1
* Ready for sid, starting transition.
* All symfiles updated to latest builds.
* Added dh_numpy call in debian/rules to depend on numpy ABI.
* Policy bumped to 3.9.3, no changes required.

Show diffs side-by-side

added added

removed removed

Lines of Context:
frmts/pcidsk/sdk/segment/cpcidskgcp2segment.cpp
 
1
/******************************************************************************
 
2
 *
 
3
 * Purpose: Implementation of access to a PCIDSK GCP2 Segment
 
4
 * 
 
5
 ******************************************************************************
 
6
 * Copyright (c) 2009
 
7
 * PCI Geomatics, 50 West Wilmot Street, Richmond Hill, Ont, Canada
 
8
 *
 
9
 * Permission is hereby granted, free of charge, to any person obtaining a
 
10
 * copy of this software and associated documentation files (the "Software"),
 
11
 * to deal in the Software without restriction, including without limitation
 
12
 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
 
13
 * and/or sell copies of the Software, and to permit persons to whom the
 
14
 * Software is furnished to do so, subject to the following conditions:
 
15
 *
 
16
 * The above copyright notice and this permission notice shall be included
 
17
 * in all copies or substantial portions of the Software.
 
18
 *
 
19
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
 
20
 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 
21
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
 
22
 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 
23
 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
 
24
 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
 
25
 * DEALINGS IN THE SOFTWARE.
 
26
 ****************************************************************************/
 
27
#include "segment/cpcidskgcp2segment.h"
 
28
 
 
29
#include "pcidsk_gcp.h"
 
30
#include "pcidsk_exception.h"
 
31
 
 
32
#include <cstring>
 
33
#include <iostream>
 
34
#include <vector>
 
35
#include <string>
 
36
 
 
37
using namespace PCIDSK;
 
38
 
 
39
struct CPCIDSKGCP2Segment::PCIDSKGCP2SegInfo
 
40
{
 
41
    std::vector<PCIDSK::GCP> gcps;
 
42
    unsigned int num_gcps;
 
43
    PCIDSKBuffer seg_data;
 
44
    
 
45
    std::string map_units;   ///< PCI mapunits string
 
46
    std::string proj_parms;  ///< Additional projection parameters
 
47
    unsigned int num_proj;
 
48
    bool changed;
 
49
};
 
50
 
 
51
CPCIDSKGCP2Segment::CPCIDSKGCP2Segment(PCIDSKFile *file, int segment, const char *segment_pointer)
 
52
    : CPCIDSKSegment(file, segment, segment_pointer), loaded_(false)
 
53
{
 
54
    pimpl_ = new PCIDSKGCP2SegInfo;
 
55
    pimpl_->gcps.clear();
 
56
    pimpl_->changed = false;
 
57
    Load();
 
58
}
 
59
 
 
60
CPCIDSKGCP2Segment::~CPCIDSKGCP2Segment()
 
61
{
 
62
    RebuildSegmentData();
 
63
    delete pimpl_;
 
64
}
 
65
 
 
66
void CPCIDSKGCP2Segment::Load()
 
67
{
 
68
    if (loaded_) {
 
69
        return;
 
70
    }
 
71
    
 
72
    // Read the the segment in. The first block has information about
 
73
    // the structure of the GCP segment (how many, the projection, etc.)
 
74
    pimpl_->seg_data.SetSize(data_size - 1024);
 
75
    ReadFromFile(pimpl_->seg_data.buffer, 0, data_size - 1024);
 
76
    
 
77
    // check for 'GCP2    ' in the first 8 bytes
 
78
    if (std::strncmp(pimpl_->seg_data.buffer, "GCP2    ", 8) != 0) {
 
79
        // Assume it's an empty segment, so we can mark loaded_ = true,
 
80
        // write it out and return
 
81
        pimpl_->changed = true;
 
82
        pimpl_->map_units = "LAT/LONG D000";
 
83
        pimpl_->proj_parms = "";
 
84
        pimpl_->num_gcps = 0;
 
85
        loaded_ = true;
 
86
        return;
 
87
    }
 
88
    
 
89
    // Check the number of blocks field's validity
 
90
    unsigned int num_blocks = pimpl_->seg_data.GetInt(8, 8);
 
91
    
 
92
    if (((data_size - 1024 - 512) / 512) != num_blocks) {
 
93
        //ThrowPCIDSKException("Calculated number of blocks (%d) does not match "
 
94
        //    "the value encoded in the GCP2 segment (%d).", ((data_size - 1024 - 512)/512), 
 
95
        //    num_blocks);
 
96
        // Something is messed up with how GDB generates these segments... nice.
 
97
    }
 
98
    
 
99
    pimpl_->num_gcps = pimpl_->seg_data.GetInt(16, 8);
 
100
    
 
101
    // Extract the map units string:
 
102
    pimpl_->map_units = std::string(pimpl_->seg_data.buffer + 24, 16);
 
103
    
 
104
    // Extract the projection parameters string
 
105
    pimpl_->proj_parms = std::string(pimpl_->seg_data.buffer + 256, 256);
 
106
    
 
107
    // Get the number of alternative projections (should be 0!)
 
108
    pimpl_->num_proj = pimpl_->seg_data.GetInt(40, 8);
 
109
    if (pimpl_->num_proj != 0) {
 
110
        ThrowPCIDSKException("There are alternative projections contained in this "
 
111
            "GCP2 segment. This functionality is not supported in libpcidsk.");
 
112
    }
 
113
    
 
114
    // Load the GCPs into the vector of PCIDSK::GCPs
 
115
    for (unsigned int i = 0; i < pimpl_->num_gcps; i++)
 
116
    {
 
117
        unsigned int offset = 512 + i * 256;
 
118
        bool is_cp = pimpl_->seg_data.buffer[offset] == 'C';
 
119
        double pixel = pimpl_->seg_data.GetDouble(offset + 6, 14);
 
120
        double line = pimpl_->seg_data.GetDouble(offset + 20, 14);
 
121
        
 
122
        double elev = pimpl_->seg_data.GetDouble(offset + 34, 12);
 
123
        double x = pimpl_->seg_data.GetDouble(offset + 48, 22);
 
124
        double y = pimpl_->seg_data.GetDouble(offset + 70, 22);
 
125
        
 
126
        PCIDSK::GCP::EElevationDatum elev_datum = pimpl_->seg_data.buffer[offset + 47] != 'M' ? 
 
127
            GCP::EEllipsoidal : GCP::EMeanSeaLevel;
 
128
        
 
129
        char elev_unit_c = pimpl_->seg_data.buffer[offset + 46];
 
130
        PCIDSK::GCP::EElevationUnit elev_unit = elev_unit_c == 'M' ? GCP::EMetres :
 
131
            elev_unit_c == 'F' ? GCP::EInternationalFeet :
 
132
            elev_unit_c == 'A' ? GCP::EAmericanFeet : GCP::EUnknown;
 
133
        
 
134
        double pix_err = pimpl_->seg_data.GetDouble(offset + 92, 10);
 
135
        double line_err = pimpl_->seg_data.GetDouble(offset + 102, 10);
 
136
        double elev_err = pimpl_->seg_data.GetDouble(offset + 112, 10);
 
137
        
 
138
        double x_err = pimpl_->seg_data.GetDouble(offset + 122, 14);
 
139
        double y_err = pimpl_->seg_data.GetDouble(offset + 136, 14);
 
140
        
 
141
        std::string gcp_id(pimpl_->seg_data.buffer + offset + 192, 64);
 
142
        
 
143
        PCIDSK::GCP gcp(x, y, elev,
 
144
                        line, pixel, gcp_id, pimpl_->map_units, 
 
145
                        pimpl_->proj_parms,
 
146
                        x_err, y_err, elev_err,
 
147
                        line_err, pix_err);
 
148
        gcp.SetElevationUnit(elev_unit);
 
149
        gcp.SetElevationDatum(elev_datum);  
 
150
        gcp.SetCheckpoint(is_cp);          
 
151
        
 
152
        pimpl_->gcps.push_back(gcp);
 
153
    } 
 
154
    
 
155
    loaded_ = true;
 
156
}
 
157
 
 
158
 // Return all GCPs in the segment
 
159
std::vector<PCIDSK::GCP> const& CPCIDSKGCP2Segment::GetGCPs(void) const
 
160
{
 
161
    return pimpl_->gcps;
 
162
}
 
163
 
 
164
// Write the given GCPs to the segment. If the segment already
 
165
// exists, it will be replaced with this one.
 
166
void CPCIDSKGCP2Segment::SetGCPs(std::vector<PCIDSK::GCP> const& gcps)
 
167
{
 
168
    pimpl_->num_gcps = gcps.size();
 
169
    pimpl_->gcps = gcps; // copy them in
 
170
    pimpl_->changed = true;
 
171
    
 
172
    RebuildSegmentData();
 
173
}
 
174
 
 
175
// Return the count of GCPs in the segment
 
176
unsigned int  CPCIDSKGCP2Segment::GetGCPCount(void) const
 
177
{
 
178
    return pimpl_->num_gcps;
 
179
}
 
180
 
 
181
void CPCIDSKGCP2Segment::RebuildSegmentData(void)
 
182
{
 
183
    if (pimpl_->changed == false) {
 
184
        return;
 
185
    }
 
186
    
 
187
    // Rebuild the segment data based on the contents of the struct
 
188
    int num_blocks = (pimpl_->num_gcps + 1) / 2;
 
189
    
 
190
    // This will have to change when we have proper projections support
 
191
 
 
192
    if (pimpl_->gcps.size() > 0)
 
193
    {
 
194
        pimpl_->gcps[0].GetMapUnits(pimpl_->map_units, 
 
195
            pimpl_->proj_parms);
 
196
    }
 
197
    
 
198
    pimpl_->seg_data.SetSize(num_blocks * 512 + 512);
 
199
    
 
200
    // Write out the first few fields
 
201
    pimpl_->seg_data.Put("GCP2    ", 0, 8);
 
202
    pimpl_->seg_data.Put(num_blocks, 8, 8);
 
203
    pimpl_->seg_data.Put((int)pimpl_->gcps.size(), 16, 8);
 
204
    pimpl_->seg_data.Put(pimpl_->map_units.c_str(), 24, 16);
 
205
    pimpl_->seg_data.Put((int)0, 40, 8);
 
206
    pimpl_->seg_data.Put(pimpl_->proj_parms.c_str(), 256, 256);
 
207
    
 
208
    // Time to write GCPs out:
 
209
    std::vector<PCIDSK::GCP>::const_iterator iter =
 
210
        pimpl_->gcps.begin();
 
211
    
 
212
    unsigned int id = 0;
 
213
    while (iter != pimpl_->gcps.end()) {
 
214
        std::size_t offset = 512 + id * 256;
 
215
        
 
216
        if ((*iter).IsCheckPoint()) {
 
217
            pimpl_->seg_data.Put("C", offset, 1);
 
218
        } else {
 
219
            pimpl_->seg_data.Put("G", offset, 1);
 
220
        }
 
221
        
 
222
        pimpl_->seg_data.Put("0", offset + 1, 5);
 
223
        
 
224
        // Start writing out the GCP values
 
225
        pimpl_->seg_data.Put((*iter).GetPixel(), offset + 6, 14, "%14.4f");
 
226
        pimpl_->seg_data.Put((*iter).GetLine(), offset + 20, 14, "%14.4f");
 
227
        pimpl_->seg_data.Put((*iter).GetZ(), offset + 34, 12, "%12.4f");
 
228
        
 
229
        GCP::EElevationUnit unit;
 
230
        GCP::EElevationDatum datum;
 
231
        (*iter).GetElevationInfo(datum, unit);
 
232
        
 
233
        char unit_c[2];
 
234
        
 
235
        switch (unit)
 
236
        {
 
237
        case GCP::EMetres:
 
238
        case GCP::EUnknown:
 
239
            unit_c[0] = 'M';
 
240
            break;
 
241
        case GCP::EAmericanFeet:
 
242
            unit_c[0] = 'A';
 
243
            break;
 
244
        case GCP::EInternationalFeet:
 
245
            unit_c[0] = 'F';
 
246
            break;
 
247
        }
 
248
        
 
249
        char datum_c[2];
 
250
        
 
251
        switch(datum)
 
252
        {
 
253
        case GCP::EEllipsoidal:
 
254
            datum_c[0] = 'E';
 
255
            break;
 
256
        case GCP::EMeanSeaLevel:
 
257
            datum_c[0] = 'M';
 
258
            break;
 
259
        }
 
260
      
 
261
        unit_c[1] = '\0';
 
262
        datum_c[1] = '\0';
 
263
        
 
264
        // Write out elevation information
 
265
        pimpl_->seg_data.Put(unit_c, offset + 46, 1);
 
266
        pimpl_->seg_data.Put(datum_c, offset + 47, 1);
 
267
        
 
268
        pimpl_->seg_data.Put((*iter).GetX(), offset + 48, 22, "%22.14e");
 
269
        pimpl_->seg_data.Put((*iter).GetY(), offset + 70, 22, "%22.14e");
 
270
        pimpl_->seg_data.Put((*iter).GetPixelErr(), offset + 92, 10, "%10.4f");
 
271
        pimpl_->seg_data.Put((*iter).GetLineErr(), offset + 102, 10, "%10.4f");
 
272
        pimpl_->seg_data.Put((*iter).GetZErr(), offset + 112, 10, "%10.4f");
 
273
        pimpl_->seg_data.Put((*iter).GetXErr(), offset + 122, 14, "%14.4e");
 
274
        pimpl_->seg_data.Put((*iter).GetYErr(), offset + 136, 14, "%14.4e");
 
275
        pimpl_->seg_data.Put((*iter).GetIDString(), offset + 192, 64, true );
 
276
        
 
277
        id++;
 
278
        iter++;
 
279
    }
 
280
    
 
281
    WriteToFile(pimpl_->seg_data.buffer, 0, pimpl_->seg_data.buffer_size);
 
282
    
 
283
    pimpl_->changed = false;
 
284
}
 
285
 
 
286
// Clear a GCP Segment
 
287
void  CPCIDSKGCP2Segment::ClearGCPs(void)
 
288
{
 
289
    pimpl_->num_gcps = 0;
 
290
    pimpl_->gcps.clear();
 
291
    pimpl_->changed = true;
 
292
    
 
293
    RebuildSegmentData();
 
294
}
 
295