source: trunk/zoo-project/zoo-services/gdal/translate/service.c @ 372

Last change on this file since 372 was 372, checked in by djay, 10 years ago

Fix for issue when using Gdal_Translate service without providing any GCP.

  • Property svn:eol-style set to native
  • Property svn:mime-type set to text/x-csrc
File size: 30.2 KB
Line 
1/******************************************************************************
2 * $Id$
3 *
4 * Project:  GDAL Utilities
5 * Purpose:  GDAL Image Translator Program
6 * Author:   Frank Warmerdam, warmerdam@pobox.com
7 *
8 ******************************************************************************
9 * Copyright (c) 1998, 2002, Frank Warmerdam
10 *
11 * Permission is hereby granted, free of charge, to any person obtaining a
12 * copy of this software and associated documentation files (the "Software"),
13 * to deal in the Software without restriction, including without limitation
14 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 * and/or sell copies of the Software, and to permit persons to whom the
16 * Software is furnished to do so, subject to the following conditions:
17 *
18 * The above copyright notice and this permission notice shall be included
19 * in all copies or substantial portions of the Software.
20 *
21 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 * DEALINGS IN THE SOFTWARE.
28 ****************************************************************************/
29
30#include "cpl_vsi.h"
31#include "cpl_conv.h"
32#include "cpl_string.h"
33#include "gdal_priv.h"
34#include "ogr_spatialref.h"
35#include "vrt/vrtdataset.h"
36
37#include "service.h"
38
39CPL_CVSID("$Id$");
40
41extern "C" {
42
43  static void AttachMetadata( GDALDatasetH, char ** );
44  static int bSubCall = FALSE;
45
46  /************************************************************************/
47  /*                          Gdal_Translate()                            */
48  /************************************************************************/
49
50#ifdef WIN32
51__declspec(dllexport)
52#endif
53  int Gdal_Translate(maps*& conf,maps*& inputs,maps*& outputs)
54  {
55
56
57    GDALDatasetH        hDataset, hOutDS;
58    int                 i;
59    int                 nRasterXSize, nRasterYSize;
60    const char          *pszSource=NULL, *pszDest=NULL, *pszFormat = "GTiff";
61    GDALDriverH         hDriver;
62    int                 *panBandList = NULL, nBandCount = 0, bDefBands = TRUE;
63    double              adfGeoTransform[6];
64    GDALDataType        eOutputType = GDT_Unknown;
65    int                 nOXSize = 0, nOYSize = 0;
66    char                *pszOXSize=NULL, *pszOYSize=NULL;
67    char                **papszCreateOptions = NULL;
68    int                 anSrcWin[4], bStrict = FALSE;
69    const char          *pszProjection;
70    int                 bScale = FALSE, bHaveScaleSrc = FALSE;
71    double              dfScaleSrcMin=0.0, dfScaleSrcMax=255.0;
72    double              dfScaleDstMin=0.0, dfScaleDstMax=255.0;
73    double              dfULX, dfULY, dfLRX, dfLRY;
74    char                **papszMetadataOptions = NULL;
75    char                *pszOutputSRS = NULL;
76    int                 bQuiet = TRUE, bGotBounds = FALSE;
77    GDALProgressFunc    pfnProgress = GDALDummyProgress;
78    int                 nGCPCount = 0;
79    GDAL_GCP            *pasGCPs = NULL;
80    int                 iSrcFileArg = -1, iDstFileArg = -1;
81    int                 bCopySubDatasets = FALSE;
82    double              adfULLR[4] = { 0,0,0,0 };
83    int                 bSetNoData = FALSE;
84    double              dfNoDataReal = 0.0;
85    int                 nRGBExpand = 0;
86
87    anSrcWin[0] = 0;
88    anSrcWin[1] = 0;
89    anSrcWin[2] = 0;
90    anSrcWin[3] = 0;
91
92    dfULX = dfULY = dfLRX = dfLRY = 0.0;
93
94    /* ----------------------------------------------------------------- */
95    /*      Register standard GDAL drivers, and process generic GDAL     */
96    /* ----------------------------------------------------------------- */
97    GDALAllRegister();
98    /* ----------------------------------------------------------------- */
99    /* Extract Format, InputDSN, OutputDSN parameters                    */
100    /* ----------------------------------------------------------------- */
101
102    map* tmpMap=NULL;
103
104    char dataPath[1024];
105    tmpMap=getMapFromMaps(conf,"main","dataPath");
106    if(tmpMap!=NULL)
107      sprintf(dataPath,"%s",tmpMap->value);
108    tmpMap=NULL;
109
110    char tempPath[1024];
111    tmpMap=getMapFromMaps(conf,"main","tmpPath");
112    if(tmpMap!=NULL){
113      sprintf(tempPath,"%s",tmpMap->value);
114    }
115    tmpMap=NULL;
116
117    tmpMap=getMapFromMaps(inputs,"Format","value");
118    if(tmpMap!=NULL){
119      pszFormat=tmpMap->value;
120    }
121    tmpMap=NULL;
122    tmpMap=getMapFromMaps(inputs,"InputDSN","value");
123    if(tmpMap!=NULL){
124      pszSource=(char*)CPLMalloc(sizeof(char)*(strlen(dataPath)+strlen(tmpMap->value)+4));
125      sprintf((char*)pszSource,"%s/%s",dataPath,tmpMap->value);
126    }
127    tmpMap=NULL;
128    tmpMap=getMapFromMaps(inputs,"OutputDSN","value");
129    if(tmpMap!=NULL){
130      pszDest=(char*)CPLMalloc(sizeof(char)*(strlen(tempPath)+strlen(tmpMap->value)+4));
131      char *ext=new char[4];
132      ext="tif";
133      if(strncasecmp(pszFormat,"AAIGRID",7)==0)
134        ext="csv";
135      else 
136        if(strncasecmp(pszFormat,"PNG",3)==0)
137          ext="png";
138        else
139          if(strncasecmp(pszFormat,"GIF",3)==0)
140            ext="gif";
141          else
142            if(strncasecmp(pszFormat,"JPEG",4)==0)
143              ext="jpg";
144      sprintf((char*)pszDest,"%s/%s.%s",tempPath,tmpMap->value,ext);
145      fprintf(stderr,"DEBUG pszDest : %s\n",pszDest);
146    }
147    tmpMap=NULL;
148    tmpMap=getMapFromMaps(inputs,"ProjWin","value");
149    if(tmpMap!=NULL){
150      char *tmp=tmpMap->value;
151      char *t=strtok(tmp,",");
152      int cnt=0;
153      while(t!=NULL){
154        switch(cnt){
155        case 0:
156          dfULX = atof(t);
157          break;
158        case 1:
159          dfULY = atof(t);
160          break;
161        case 2:
162          dfLRX = atof(t);
163          break;
164        case 3:
165          dfLRY = atof(t);
166          break;
167        }
168        t=strtok(NULL,",");
169        cnt++;
170      }
171    }
172
173    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
174    tmpMap=NULL;
175    tmpMap=getMapFromMaps(inputs,"GCP","value");
176    if(tmpMap!=NULL && strncasecmp(tmpMap->value,"NULL",4)!=0){
177      map* length=getMapFromMaps(inputs,"GCP","length");
178      int len=0;
179      if(length){
180        len=atoi(length->value);
181        int i;
182        maps* currentMaps=getMaps(inputs,"GCP");
183        for(i=0;i<len;i++){
184          char* endptr = NULL;
185          /* -gcp pixel line easting northing [elev] */
186         
187          nGCPCount++;
188          pasGCPs = (GDAL_GCP *) 
189            realloc( pasGCPs, sizeof(GDAL_GCP) * nGCPCount );
190          GDALInitGCPs( 1, pasGCPs + nGCPCount - 1 );
191
192          map* currentMap=getMapArray(currentMaps->content,"value",i);
193
194          char* tmpV=strdup(currentMap->value);
195          char *res=strtok(tmpV,",");
196          int j=0;
197          while(res!=NULL){
198            switch(j){
199            case 0:
200              pasGCPs[nGCPCount-1].dfGCPPixel = CPLAtofM(res);
201              break;
202            case 1:
203              pasGCPs[nGCPCount-1].dfGCPLine = CPLAtofM(res);
204              break;
205            case 2:
206              pasGCPs[nGCPCount-1].dfGCPX = CPLAtofM(res);
207              break;
208            case 3:
209              pasGCPs[nGCPCount-1].dfGCPY = CPLAtofM(res);
210              break;
211            case 4:
212              if(res!=NULL && (strtod(res, &endptr) != 0.0 || res[0] == '0'))
213                if (endptr && *endptr == 0)
214                  pasGCPs[nGCPCount-1].dfGCPZ = CPLAtofM(res);
215              break;
216            }
217            res=strtok(NULL,",");
218            j++;
219          }
220        }
221      }else{
222        char* endptr = NULL;
223        /* -gcp pixel line easting northing [elev] */
224       
225        nGCPCount++;
226        pasGCPs = (GDAL_GCP *) 
227          realloc( pasGCPs, sizeof(GDAL_GCP) * nGCPCount );
228        GDALInitGCPs( 1, pasGCPs + nGCPCount - 1 );
229       
230        char* tmpV=strdup(tmpMap->value);
231        char *res=strtok(tmpV,",");
232        int j=0;
233        while(res!=NULL){
234          switch(j){
235          case 0:
236            pasGCPs[nGCPCount-1].dfGCPPixel = CPLAtofM(res);
237            break;
238          case 1:
239            pasGCPs[nGCPCount-1].dfGCPLine = CPLAtofM(res);
240            break;
241          case 2:
242            pasGCPs[nGCPCount-1].dfGCPX = CPLAtofM(res);
243            break;
244          case 3:
245            pasGCPs[nGCPCount-1].dfGCPY = CPLAtofM(res);
246            break;
247          case 4:
248            if(res!=NULL && (CPLStrtod(res, &endptr) != 0.0 || res[0] == '0'))
249              if (endptr && *endptr == 0)
250                pasGCPs[nGCPCount-1].dfGCPZ = CPLAtofM(res);
251            break;
252          }         
253          res=strtok(NULL,",");
254          j++;
255        }
256      }
257    }
258
259    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
260    tmpMap=NULL;
261    tmpMap=getMapFromMaps(inputs,"SRS","value");
262    if(tmpMap!=NULL){
263      OGRSpatialReference oOutputSRS;
264      if( oOutputSRS.SetFromUserInput( tmpMap->value ) != OGRERR_NONE )
265        {
266          char *msg=(char*)CPLMalloc(100*sizeof(char));
267          sprintf( msg, "Failed to process SRS definition: %s\n", 
268                   tmpMap->value );
269          setMapInMaps(conf,"lenv","message",msg);
270          /**
271           * Avoiding GDALDestroyDriverManager() call
272           */
273          return SERVICE_FAILED;
274        }
275      oOutputSRS.exportToWkt( &pszOutputSRS );
276    }
277    tmpMap=NULL;
278    tmpMap=getMapFromMaps(inputs,"Type","value");
279    if(tmpMap!=NULL){
280      int       iType;
281     
282      for( iType = 1; iType < GDT_TypeCount; iType++ )
283        {
284          if( GDALGetDataTypeName((GDALDataType)iType) != NULL
285              && EQUAL(GDALGetDataTypeName((GDALDataType)iType),
286                       tmpMap->value) )
287            {
288              eOutputType = (GDALDataType) iType;
289            }
290        }
291     
292      if( eOutputType == GDT_Unknown )
293        {
294          printf( "Unknown output pixel type: %s\n", tmpMap->value );
295          /**
296           * Avoiding GDALDestroyDriverManager() call
297           */
298          exit( 2 );
299        }
300    }
301
302    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
303    if( pszDest == NULL ){
304        fprintf(stderr,"exit line 416");
305        fflush(stderr);
306        /**
307         * Avoiding GDALDestroyDriverManager() call
308         */
309        exit( 10 );
310      }
311
312    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
313    if ( strcmp(pszSource, pszDest) == 0)
314      {
315        fprintf(stderr, "Source and destination datasets must be different.\n");
316        fflush(stderr);
317        /**
318         * Avoiding GDALDestroyDriverManager() call
319         */
320        exit( 1 );
321      }
322
323    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
324    /* ----------------------------------------------------------------- */
325    /*      Attempt to open source file.                                 */
326    /* ----------------------------------------------------------------- */
327
328    hDataset = GDALOpenShared( pszSource, GA_ReadOnly );
329    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
330   
331    if( hDataset == NULL ){
332      char *msg=(char*) CPLMalloc(1024*sizeof(char));
333      sprintf( msg,
334               "GDALOpen failed - %d\n%s\n",
335               CPLGetLastErrorNo(), CPLGetLastErrorMsg() );
336      setMapInMaps(conf,"lenv","message",msg);
337      return SERVICE_FAILED;
338    }
339
340    fprintf(stderr,"DEBUG pszDest : %s %d\n",pszDest,__LINE__);
341    /* ----------------------------------------------------------------- */
342    /*      Handle subdatasets.                                          */
343    /* ----------------------------------------------------------------- */
344    if( !bCopySubDatasets
345        && CSLCount(GDALGetMetadata( hDataset, "SUBDATASETS" )) > 0 
346        && GDALGetRasterCount(hDataset) == 0 )
347      {
348        char *msg=(char*) CPLMalloc(1024*sizeof(char));
349        sprintf( msg,
350                 "Input file contains subdatasets. Please, select one of them for reading.\n" );
351        setMapInMaps(conf,"lenv","message",msg);
352        GDALClose( hDataset );
353        /**
354         * Avoiding GDALDestroyDriverManager() call
355         */
356        exit( 1 );
357      }
358
359    if( CSLCount(GDALGetMetadata( hDataset, "SUBDATASETS" )) > 0 
360        && bCopySubDatasets )
361      {
362        char **papszSubdatasets = GDALGetMetadata(hDataset,"SUBDATASETS");
363        char *pszSubDest = (char *) CPLMalloc(strlen(pszDest)+32);
364        int i;
365        int bOldSubCall = bSubCall;
366       
367        //argv[iDstFileArg] = pszSubDest;
368        bSubCall = TRUE;
369        for( i = 0; papszSubdatasets[i] != NULL; i += 2 )
370          {
371            //argv[iSrcFileArg] = strstr(papszSubdatasets[i],"=")+1;
372            sprintf( pszSubDest, "%s%d", pszDest, i/2 + 1 );
373            /*if( ProxyMain( argc, argv ) != 0 )
374              break;*/
375          }
376       
377        bSubCall = bOldSubCall;
378        free( pszSubDest );
379
380        GDALClose( hDataset );
381
382        if( !bSubCall )
383          {
384            GDALDumpOpenDatasets( stderr );
385            fflush(stderr);
386            /**
387             * Avoiding GDALDestroyDriverManager() call
388             */
389          }
390        return 1;
391      }
392
393    /* ----------------------------------------------------------------- */
394    /*      Collect some information from the source file.               */
395    /* ----------------------------------------------------------------- */
396    nRasterXSize = GDALGetRasterXSize( hDataset );
397    nRasterYSize = GDALGetRasterYSize( hDataset );
398
399    if( !bQuiet )
400      fprintf( stderr, "Input file size is %d, %d\n", nRasterXSize, nRasterYSize );
401
402    if( anSrcWin[2] == 0 && anSrcWin[3] == 0 ){
403        anSrcWin[2] = nRasterXSize;
404        anSrcWin[3] = nRasterYSize;
405      }
406
407    /* ----------------------------------------------------------------- */
408    /*  Build band list to translate                                     */
409    /* ----------------------------------------------------------------- */
410    if( nBandCount == 0 ){
411        nBandCount = GDALGetRasterCount( hDataset );
412        if( nBandCount == 0 ){
413            fprintf( stderr, "Input file has no bands, and so cannot be translated.\n" );
414            fflush(stderr);
415            /**
416             * Avoiding GDALDestroyDriverManager() call
417             */
418            exit(1 );
419          }
420
421        panBandList = (int *) CPLMalloc(sizeof(int)*nBandCount);
422        for( i = 0; i < nBandCount; i++ )
423          panBandList[i] = i+1;
424      }
425    else
426      {
427        for( i = 0; i < nBandCount; i++ )
428          {
429            if( panBandList[i] < 1 || panBandList[i] > GDALGetRasterCount(hDataset) )
430              {
431                fprintf( stderr, 
432                         "Band %d requested, but only bands 1 to %d available.\n",
433                         panBandList[i], GDALGetRasterCount(hDataset) );
434                fflush(stderr);
435                /**
436                 * Avoiding GDALDestroyDriverManager() call
437                 */
438                exit( 2 );
439              }
440          }
441
442        if( nBandCount != GDALGetRasterCount( hDataset ) )
443          bDefBands = FALSE;
444      }
445
446    /* ----------------------------------------------------------------- */
447    /*   Compute the source window from the projected source window      */
448    /*   if the projected coordinates were provided.  Note that the      */
449    /*   projected coordinates are in ulx, uly, lrx, lry format,         */
450    /*   while the anSrcWin is xoff, yoff, xsize, ysize with the         */
451    /*   xoff,yoff being the ulx, uly in pixel/line.                     */
452    /* ----------------------------------------------------------------- */
453    if( dfULX != 0.0 || dfULY != 0.0 
454        || dfLRX != 0.0 || dfLRY != 0.0 )
455      {
456        double  adfGeoTransform[6];
457
458        GDALGetGeoTransform( hDataset, adfGeoTransform );
459
460        if( adfGeoTransform[2] != 0.0 || adfGeoTransform[4] != 0.0 ){
461            fprintf( stderr, 
462                     "The -projwin option was used, but the geotransform is\n"
463                     "rotated.  This configuration is not supported.\n" );
464            GDALClose( hDataset );
465            free( panBandList );
466            fflush(stderr);
467            /**
468             * Avoiding GDALDestroyDriverManager() call
469             */
470            exit( 1 );
471          }
472
473        anSrcWin[0] = (int) 
474          ((dfULX - adfGeoTransform[0]) / adfGeoTransform[1] + 0.001);
475        anSrcWin[1] = (int) 
476          ((dfULY - adfGeoTransform[3]) / adfGeoTransform[5] + 0.001);
477
478        anSrcWin[2] = (int) ((dfLRX - dfULX) / adfGeoTransform[1] + 0.5);
479        anSrcWin[3] = (int) ((dfLRY - dfULY) / adfGeoTransform[5] + 0.5);
480
481        if( !bQuiet )
482          fprintf( stdout, 
483                   "Computed -srcwin %d %d %d %d from projected window.\n",
484                   anSrcWin[0], 
485                   anSrcWin[1], 
486                   anSrcWin[2], 
487                   anSrcWin[3] );
488       
489        if( anSrcWin[0] < 0 || anSrcWin[1] < 0 
490            || anSrcWin[0] + anSrcWin[2] > GDALGetRasterXSize(hDataset) 
491            || anSrcWin[1] + anSrcWin[3] > GDALGetRasterYSize(hDataset) )
492          {
493            fprintf( stderr, 
494                     "Computed -srcwin falls outside raster size of %dx%d.\n",
495                     GDALGetRasterXSize(hDataset), 
496                     GDALGetRasterYSize(hDataset) );
497            exit( 1 );
498          }
499      }
500
501    /* ----------------------------------------------------------------- */
502    /*      Verify source window.                                        */
503    /* ----------------------------------------------------------------- */
504    if( anSrcWin[0] < 0 || anSrcWin[1] < 0 
505        || anSrcWin[2] <= 0 || anSrcWin[3] <= 0
506        || anSrcWin[0] + anSrcWin[2] > GDALGetRasterXSize(hDataset) 
507        || anSrcWin[1] + anSrcWin[3] > GDALGetRasterYSize(hDataset) )
508      {
509        fprintf( stderr, 
510                 "-srcwin %d %d %d %d falls outside raster size of %dx%d\n"
511                 "or is otherwise illegal.\n",
512                 anSrcWin[0],
513                 anSrcWin[1],
514                 anSrcWin[2],
515                 anSrcWin[3],
516                 GDALGetRasterXSize(hDataset), 
517                 GDALGetRasterYSize(hDataset) );
518        exit( 1 );
519      }
520
521    /* ----------------------------------------------------------------- */
522    /*      Find the output driver.                                      */
523    /* ----------------------------------------------------------------- */
524    hDriver = GDALGetDriverByName( pszFormat );
525    if( hDriver == NULL )
526      {
527        int     iDr;
528       
529        char* msg=(char*) CPLMalloc(4096*sizeof(char));
530        sprintf(msg,"Output driver `%s' not recognised.\nThe following format drivers are configured and support output:\n",pszFormat);
531        for( iDr = 0; iDr < GDALGetDriverCount(); iDr++ )
532          {
533            GDALDriverH hDriver = GDALGetDriver(iDr);
534
535            if( GDALGetMetadataItem( hDriver, GDAL_DCAP_CREATE, NULL ) != NULL
536                || GDALGetMetadataItem( hDriver, GDAL_DCAP_CREATECOPY,
537                                        NULL ) != NULL )
538              {
539                fprintf(stderr,msg);
540                char *tmp=strdup(msg);
541                sprintf( msg,"%s  %s: %s\n",tmp,
542                        GDALGetDriverShortName( hDriver  ),
543                        GDALGetDriverLongName( hDriver ) );
544                free(tmp);
545              }
546          }
547        setMapInMaps(conf,"lenv","message",msg);
548        GDALClose( hDataset );
549        free( panBandList );
550        fflush(stderr);
551        /**
552         * Avoiding GDALDestroyDriverManager() call
553         */
554        CSLDestroy( papszCreateOptions );
555        return 4;
556      }
557
558    /* ----------------------------------------------------------------- */
559    /*   The short form is to CreateCopy().  We use this if the input    */
560    /*   matches the whole dataset.  Eventually we should rewrite        */
561    /*   this entire program to use virtual datasets to construct a      */
562    /*   virtual input source to copy from.                              */
563    /* ----------------------------------------------------------------- */
564    int bSpatialArrangementPreserved = (
565           anSrcWin[0] == 0 && anSrcWin[1] == 0
566        && anSrcWin[2] == GDALGetRasterXSize(hDataset)
567        && anSrcWin[3] == GDALGetRasterYSize(hDataset)
568        && pszOXSize == NULL && pszOYSize == NULL );
569
570    if( eOutputType == GDT_Unknown
571        && !bScale && CSLCount(papszMetadataOptions) == 0 && bDefBands
572        && bSpatialArrangementPreserved
573        && nGCPCount == 0 && !bGotBounds
574        && pszOutputSRS == NULL && !bSetNoData
575        && nRGBExpand == 0)
576      {
577       
578        hOutDS = GDALCreateCopy( hDriver, pszDest, hDataset, 
579                                 FALSE, (char**)0, 
580                                 pfnProgress, NULL );
581
582        if( hOutDS != NULL )
583          GDALClose( hOutDS );
584       
585        GDALClose( hDataset );
586
587        free( panBandList );
588
589        if( !bSubCall )
590          {
591            GDALDumpOpenDatasets( stderr );
592            /**
593             * Avoiding GDALDestroyDriverManager() call
594             */
595          }
596        fprintf(stderr,"==%s %s %s %d==\n",pszFormat,pszSource,pszDest,__LINE__);
597        fflush(stderr);
598
599        CSLDestroy( papszCreateOptions );
600        fprintf(stderr,"==%s %s %s %d==\n",pszFormat,pszSource,pszDest,__LINE__);
601        fflush(stderr);
602        //outputs=(maps*)CPLMalloc(sizeof(maps*));
603        //outputs->name="OutputedPolygon";
604        //outputs->content=createMap("value",(char*)pszDest);
605        //outputs->next=NULL;
606        dumpMaps(outputs);
607        setMapInMaps(outputs,"Result","value",(char*)pszDest);
608        return SERVICE_SUCCEEDED;
609      }
610
611    /* ----------------------------------------------------------------- */
612    /*      Establish some parameters.                                   */
613    /* ----------------------------------------------------------------- */
614    if( pszOXSize == NULL )
615      {
616        nOXSize = anSrcWin[2];
617        nOYSize = anSrcWin[3];
618      }
619    else
620      {
621        nOXSize = (int) ((pszOXSize[strlen(pszOXSize)-1]=='%' 
622                          ? atof(pszOXSize)/100*anSrcWin[2] : atoi(pszOXSize)));
623        nOYSize = (int) ((pszOYSize[strlen(pszOYSize)-1]=='%' 
624                          ? atof(pszOYSize)/100*anSrcWin[3] : atoi(pszOYSize)));
625      }
626
627    /* ================================================================= */
628    /*      Create a virtual dataset.                                    */
629    /* ================================================================= */
630    VRTDataset *poVDS;
631       
632    /* ----------------------------------------------------------------- */
633    /*      Make a virtual clone.                                        */
634    /* ----------------------------------------------------------------- */
635    poVDS = (VRTDataset *) VRTCreate( nOXSize, nOYSize );
636
637    if( nGCPCount == 0 )
638      {
639        if( pszOutputSRS != NULL )
640          {
641            poVDS->SetProjection( pszOutputSRS );
642          }
643        else
644          {
645            pszProjection = GDALGetProjectionRef( hDataset );
646            if( pszProjection != NULL && strlen(pszProjection) > 0 )
647              poVDS->SetProjection( pszProjection );
648          }
649      }
650
651    if( bGotBounds )
652      {
653        adfGeoTransform[0] = adfULLR[0];
654        adfGeoTransform[1] = (adfULLR[2] - adfULLR[0]) / nOXSize;
655        adfGeoTransform[2] = 0.0;
656        adfGeoTransform[3] = adfULLR[1];
657        adfGeoTransform[4] = 0.0;
658        adfGeoTransform[5] = (adfULLR[3] - adfULLR[1]) / nOYSize;
659
660        poVDS->SetGeoTransform( adfGeoTransform );
661      }
662
663    else if( GDALGetGeoTransform( hDataset, adfGeoTransform ) == CE_None
664             && nGCPCount == 0 )
665      {
666        adfGeoTransform[0] += anSrcWin[0] * adfGeoTransform[1]
667          + anSrcWin[1] * adfGeoTransform[2];
668        adfGeoTransform[3] += anSrcWin[0] * adfGeoTransform[4]
669          + anSrcWin[1] * adfGeoTransform[5];
670       
671        adfGeoTransform[1] *= anSrcWin[2] / (double) nOXSize;
672        adfGeoTransform[2] *= anSrcWin[3] / (double) nOYSize;
673        adfGeoTransform[4] *= anSrcWin[2] / (double) nOXSize;
674        adfGeoTransform[5] *= anSrcWin[3] / (double) nOYSize;
675       
676        poVDS->SetGeoTransform( adfGeoTransform );
677      }
678
679    if( nGCPCount != 0 )
680      {
681        const char *pszGCPProjection = pszOutputSRS;
682
683        if( pszGCPProjection == NULL )
684          pszGCPProjection = GDALGetGCPProjection( hDataset );
685        if( pszGCPProjection == NULL )
686          pszGCPProjection = "";
687
688        poVDS->SetGCPs( nGCPCount, pasGCPs, pszGCPProjection );
689
690        GDALDeinitGCPs( nGCPCount, pasGCPs );
691        free( pasGCPs );
692      }
693
694    else if( GDALGetGCPCount( hDataset ) > 0 )
695      {
696        GDAL_GCP *pasGCPs;
697        int       nGCPs = GDALGetGCPCount( hDataset );
698
699        pasGCPs = GDALDuplicateGCPs( nGCPs, GDALGetGCPs( hDataset ) );
700
701        for( i = 0; i < nGCPs; i++ )
702          {
703            pasGCPs[i].dfGCPPixel -= anSrcWin[0];
704            pasGCPs[i].dfGCPLine  -= anSrcWin[1];
705            pasGCPs[i].dfGCPPixel *= (nOXSize / (double) anSrcWin[2] );
706            pasGCPs[i].dfGCPLine  *= (nOYSize / (double) anSrcWin[3] );
707          }
708           
709        poVDS->SetGCPs( nGCPs, pasGCPs,
710                        GDALGetGCPProjection( hDataset ) );
711
712        GDALDeinitGCPs( nGCPs, pasGCPs );
713        free( pasGCPs );
714      }
715
716    /* ----------------------------------------------------------------- */
717    /*      Transfer generally applicable metadata.                      */
718    /* ----------------------------------------------------------------- */
719    poVDS->SetMetadata( ((GDALDataset*)hDataset)->GetMetadata() );
720    AttachMetadata( (GDALDatasetH) poVDS, papszMetadataOptions );
721
722    /* ----------------------------------------------------------------- */
723    /*      Transfer metadata that remains valid if the spatial          */
724    /*      arrangement of the data is unaltered.                        */
725    /* ----------------------------------------------------------------- */
726    if( anSrcWin[0] == 0 && anSrcWin[1] == 0 
727        && anSrcWin[2] == GDALGetRasterXSize(hDataset)
728        && anSrcWin[3] == GDALGetRasterYSize(hDataset) 
729        && pszOXSize == NULL && pszOYSize == NULL )
730      {
731        char **papszMD;
732
733        papszMD = ((GDALDataset*)hDataset)->GetMetadata("RPC");
734        if( papszMD != NULL )
735          poVDS->SetMetadata( papszMD, "RPC" );
736      }
737
738    if (nRGBExpand != 0)
739      nBandCount += nRGBExpand - 1;
740
741    /* ================================================================= */
742    /*      Process all bands.                                           */
743    /* ================================================================= */
744    for( i = 0; i < nBandCount; i++ )
745      {
746        VRTSourcedRasterBand   *poVRTBand;
747        GDALRasterBand  *poSrcBand;
748        GDALDataType    eBandType;
749
750        if (nRGBExpand != 0 && i < nRGBExpand)
751          {
752            poSrcBand = ((GDALDataset *) 
753                         hDataset)->GetRasterBand(panBandList[0]);
754            if (poSrcBand->GetColorTable() == NULL)
755              {
756                fprintf(stderr, "Error : band %d has no color table\n", panBandList[0]);
757                GDALClose( hDataset );
758                free( panBandList );
759                fflush(stderr);
760                /**
761                 * Avoiding GDALDestroyDriverManager() call
762                 */
763                CSLDestroy( papszCreateOptions );
764                exit( 1 );
765              }
766          }
767        else
768          poSrcBand = ((GDALDataset *) 
769                       hDataset)->GetRasterBand(panBandList[i]);
770
771        /* ------------------------------------------------------------ */
772        /*      Select output data type to match source.                */
773        /* ------------------------------------------------------------ */
774        if( eOutputType == GDT_Unknown )
775          eBandType = poSrcBand->GetRasterDataType();
776        else
777          eBandType = eOutputType;
778
779        /* ------------------------------------------------------------ */
780        /*      Create this band.                                       */
781        /* ------------------------------------------------------------ */
782        poVDS->AddBand( eBandType, NULL );
783        poVRTBand = (VRTSourcedRasterBand *) poVDS->GetRasterBand( i+1 );
784           
785        /* ------------------------------------------------------------ */
786        /*      Do we need to collect scaling information?              */
787        /* ------------------------------------------------------------ */
788        double dfScale=1.0, dfOffset=0.0;
789
790        if( bScale && !bHaveScaleSrc )
791          {
792            double      adfCMinMax[2];
793            GDALComputeRasterMinMax( poSrcBand, TRUE, adfCMinMax );
794            dfScaleSrcMin = adfCMinMax[0];
795            dfScaleSrcMax = adfCMinMax[1];
796          }
797
798        if( bScale )
799          {
800            if( dfScaleSrcMax == dfScaleSrcMin )
801              dfScaleSrcMax += 0.1;
802            if( dfScaleDstMax == dfScaleDstMin )
803              dfScaleDstMax += 0.1;
804
805            dfScale = (dfScaleDstMax - dfScaleDstMin) 
806              / (dfScaleSrcMax - dfScaleSrcMin);
807            dfOffset = -1 * dfScaleSrcMin * dfScale + dfScaleDstMin;
808          }
809
810        /* ------------------------------------------------------------ */
811        /*      Create a simple or complex data source depending on the */
812        /*      translation type required.                              */
813        /* ------------------------------------------------------------ */
814        if( bScale || (nRGBExpand != 0 && i < nRGBExpand) )
815          {
816            poVRTBand->AddComplexSource( poSrcBand,
817                                         anSrcWin[0], anSrcWin[1], 
818                                         anSrcWin[2], anSrcWin[3], 
819                                         0, 0, nOXSize, nOYSize,
820                                         dfOffset, dfScale,
821                                         VRT_NODATA_UNSET,
822                                         (nRGBExpand != 0 && i < nRGBExpand) ? i + 1 : 0 );
823          }
824        else
825          poVRTBand->AddSimpleSource( poSrcBand,
826                                      anSrcWin[0], anSrcWin[1], 
827                                      anSrcWin[2], anSrcWin[3], 
828                                      0, 0, nOXSize, nOYSize );
829
830        /* In case of color table translate, we only set the color interpretation */
831        /* other info copied by CopyCommonInfoFrom are not relevant in RGB expansion */
832        if (nRGBExpand != 0 && i < nRGBExpand)
833          {
834            poVRTBand->SetColorInterpretation( (GDALColorInterp) (GCI_RedBand + i) );
835          }
836        else
837          {
838            /* --------------------------------------------------------- */
839            /*      copy over some other information of interest.        */
840            /* --------------------------------------------------------- */
841            poVRTBand->CopyCommonInfoFrom( poSrcBand );
842          }
843
844        /* ------------------------------------------------------------- */
845        /*      Set a forcable nodata value?                             */
846        /* ------------------------------------------------------------- */
847        if( bSetNoData )
848          poVRTBand->SetNoDataValue( dfNoDataReal );
849      }
850
851    /* ----------------------------------------------------------------- */
852    /*      Write to the output file using CopyCreate().                 */
853    /* ----------------------------------------------------------------- */
854    hOutDS = GDALCreateCopy( hDriver, pszDest, (GDALDatasetH) poVDS,
855                             bStrict, papszCreateOptions, 
856                             pfnProgress, NULL );
857
858    if( hOutDS != NULL )
859      {
860        GDALClose( hOutDS );
861      }
862   
863    GDALClose( (GDALDatasetH) poVDS );
864         
865    GDALClose( hDataset );
866   
867    free( panBandList );
868   
869    free( pszOutputSRS );
870
871    if( !bSubCall )
872      {
873        GDALDumpOpenDatasets( stderr );
874        fflush(stderr);
875        /**
876         * Avoiding GDALDestroyDriverManager() call
877         */
878      }
879
880    CSLDestroy( papszCreateOptions );
881   
882    setMapInMaps(outputs,"Result","value",(char*)pszDest);
883
884    return SERVICE_SUCCEEDED;
885  }
886
887
888  /************************************************************************/
889  /*                           AttachMetadata()                           */
890  /************************************************************************/
891
892  static void AttachMetadata( GDALDatasetH hDS, char **papszMetadataOptions )
893
894  {
895    int nCount = CSLCount(papszMetadataOptions);
896    int i;
897
898    for( i = 0; i < nCount; i++ )
899      {
900        char    *pszKey = NULL;
901        const char *pszValue;
902       
903        pszValue = CPLParseNameValue( papszMetadataOptions[i], &pszKey );
904        GDALSetMetadataItem(hDS,pszKey,pszValue,NULL);
905        free( pszKey );
906      }
907
908    CSLDestroy( papszMetadataOptions );
909  }
910
911}
Note: See TracBrowser for help on using the repository browser.

Search

ZOO Sponsors

http://www.zoo-project.org/trac/chrome/site/img/geolabs-logo.pnghttp://www.zoo-project.org/trac/chrome/site/img/neogeo-logo.png http://www.zoo-project.org/trac/chrome/site/img/apptech-logo.png http://www.zoo-project.org/trac/chrome/site/img/3liz-logo.png http://www.zoo-project.org/trac/chrome/site/img/gateway-logo.png

Become a sponsor !

Knowledge partners

http://www.zoo-project.org/trac/chrome/site/img/ocu-logo.png http://www.zoo-project.org/trac/chrome/site/img/gucas-logo.png http://www.zoo-project.org/trac/chrome/site/img/polimi-logo.png http://www.zoo-project.org/trac/chrome/site/img/fem-logo.png http://www.zoo-project.org/trac/chrome/site/img/supsi-logo.png http://www.zoo-project.org/trac/chrome/site/img/cumtb-logo.png

Become a knowledge partner

Related links

http://zoo-project.org/img/ogclogo.png http://zoo-project.org/img/osgeologo.png