LCOV - code coverage report
Current view: directory - frmts/vrt - vrtwarped.cpp (source / functions) Found Hit Coverage
Test: gdal_filtered.info Lines: 387 306 79.1 %
Date: 2012-12-26 Functions: 34 24 70.6 %

       1                 : /******************************************************************************
       2                 :  * $Id: vrtwarped.cpp 24190 2012-04-01 20:50:44Z rouault $
       3                 :  *
       4                 :  * Project:  Virtual GDAL Datasets
       5                 :  * Purpose:  Implementation of VRTWarpedRasterBand *and VRTWarpedDataset.
       6                 :  * Author:   Frank Warmerdam <warmerdam@pobox.com>
       7                 :  *
       8                 :  ******************************************************************************
       9                 :  * Copyright (c) 2004, Frank Warmerdam <warmerdam@pobox.com>
      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 "vrtdataset.h"
      31                 : #include "cpl_minixml.h"
      32                 : #include "cpl_string.h"
      33                 : #include "gdalwarper.h"
      34                 : #include "gdal_alg_priv.h"
      35                 : #include <cassert>
      36                 : 
      37                 : CPL_CVSID("$Id: vrtwarped.cpp 24190 2012-04-01 20:50:44Z rouault $");
      38                 : 
      39                 : /************************************************************************/
      40                 : /*                      GDALAutoCreateWarpedVRT()                       */
      41                 : /************************************************************************/
      42                 : 
      43                 : /**
      44                 :  * Create virtual warped dataset automatically.
      45                 :  *
      46                 :  * This function will create a warped virtual file representing the 
      47                 :  * input image warped into the target coordinate system.  A GenImgProj
      48                 :  * transformation is created to accomplish any required GCP/Geotransform
      49                 :  * warp and reprojection to the target coordinate system.  The output virtual
      50                 :  * dataset will be "northup" in the target coordinate system.   The
      51                 :  * GDALSuggestedWarpOutput() function is used to determine the bounds and
      52                 :  * resolution of the output virtual file which should be large enough to 
      53                 :  * include all the input image 
      54                 :  *
      55                 :  * Note that the constructed GDALDatasetH will acquire one or more references 
      56                 :  * to the passed in hSrcDS.  Reference counting semantics on the source 
      57                 :  * dataset should be honoured.  That is, don't just GDALClose() it unless it 
      58                 :  * was opened with GDALOpenShared(). 
      59                 :  *
      60                 :  * The returned dataset will have no associated filename for itself.  If you
      61                 :  * want to write the virtual dataset description to a file, use the
      62                 :  * GDALSetDescription() function (or SetDescription() method) on the dataset
      63                 :  * to assign a filename before it is closed.  
      64                 :  *
      65                 :  * @param hSrcDS The source dataset. 
      66                 :  *
      67                 :  * @param pszSrcWKT The coordinate system of the source image.  If NULL, it 
      68                 :  * will be read from the source image. 
      69                 :  *
      70                 :  * @param pszDstWKT The coordinate system to convert to.  If NULL no change 
      71                 :  * of coordinate system will take place.  
      72                 :  *
      73                 :  * @param eResampleAlg One of GRA_NearestNeighbour, GRA_Bilinear, GRA_Cubic or 
      74                 :  * GRA_CubicSpline.  Controls the sampling method used. 
      75                 :  *
      76                 :  * @param dfMaxError Maximum error measured in input pixels that is allowed in 
      77                 :  * approximating the transformation (0.0 for exact calculations).
      78                 :  *
      79                 :  * @param psOptionsIn Additional warp options, normally NULL.
      80                 :  *
      81                 :  * @return NULL on failure, or a new virtual dataset handle on success.
      82                 :  */
      83                 : 
      84                 : GDALDatasetH CPL_STDCALL 
      85              14 : GDALAutoCreateWarpedVRT( GDALDatasetH hSrcDS, 
      86                 :                          const char *pszSrcWKT,
      87                 :                          const char *pszDstWKT,
      88                 :                          GDALResampleAlg eResampleAlg, 
      89                 :                          double dfMaxError, 
      90                 :                          const GDALWarpOptions *psOptionsIn )
      91                 :     
      92                 : {
      93                 :     GDALWarpOptions *psWO;
      94                 :     int i;
      95                 : 
      96              14 :     VALIDATE_POINTER1( hSrcDS, "GDALAutoCreateWarpedVRT", NULL );
      97                 : 
      98                 : /* -------------------------------------------------------------------- */
      99                 : /*      Populate the warp options.                                      */
     100                 : /* -------------------------------------------------------------------- */
     101              14 :     if( psOptionsIn != NULL )
     102               0 :         psWO = GDALCloneWarpOptions( psOptionsIn );
     103                 :     else
     104              14 :         psWO = GDALCreateWarpOptions();
     105                 : 
     106              14 :     psWO->eResampleAlg = eResampleAlg;
     107                 : 
     108              14 :     psWO->hSrcDS = hSrcDS;
     109                 : 
     110              14 :     psWO->nBandCount = GDALGetRasterCount( hSrcDS );
     111              14 :     psWO->panSrcBands = (int *) CPLMalloc(sizeof(int) * psWO->nBandCount);
     112              14 :     psWO->panDstBands = (int *) CPLMalloc(sizeof(int) * psWO->nBandCount);
     113                 : 
     114              30 :     for( i = 0; i < psWO->nBandCount; i++ )
     115                 :     {
     116              16 :         psWO->panSrcBands[i] = i+1;
     117              16 :         psWO->panDstBands[i] = i+1;
     118                 :     }
     119                 : 
     120                 :     /* TODO: should fill in no data where available */
     121                 : 
     122                 : /* -------------------------------------------------------------------- */
     123                 : /*      Create the transformer.                                         */
     124                 : /* -------------------------------------------------------------------- */
     125              14 :     psWO->pfnTransformer = GDALGenImgProjTransform;
     126                 :     psWO->pTransformerArg = 
     127                 :         GDALCreateGenImgProjTransformer( psWO->hSrcDS, pszSrcWKT, 
     128                 :                                          NULL, pszDstWKT,
     129              14 :                                          TRUE, 1.0, 0 );
     130                 : 
     131              14 :     if( psWO->pTransformerArg == NULL )
     132                 :     {
     133               0 :         GDALDestroyWarpOptions( psWO );
     134               0 :         return NULL;
     135                 :     }
     136                 : 
     137                 : /* -------------------------------------------------------------------- */
     138                 : /*      Figure out the desired output bounds and resolution.            */
     139                 : /* -------------------------------------------------------------------- */
     140                 :     double adfDstGeoTransform[6];
     141                 :     int    nDstPixels, nDstLines;
     142                 :     CPLErr eErr;
     143                 : 
     144                 :     eErr = 
     145                 :         GDALSuggestedWarpOutput( hSrcDS, psWO->pfnTransformer, 
     146                 :                                  psWO->pTransformerArg, 
     147              14 :                                  adfDstGeoTransform, &nDstPixels, &nDstLines );
     148                 : 
     149                 : /* -------------------------------------------------------------------- */
     150                 : /*      Update the transformer to include an output geotransform        */
     151                 : /*      back to pixel/line coordinates.                                 */
     152                 : /*                                                                      */
     153                 : /* -------------------------------------------------------------------- */
     154                 :     GDALSetGenImgProjTransformerDstGeoTransform( 
     155              14 :         psWO->pTransformerArg, adfDstGeoTransform );
     156                 : 
     157                 : /* -------------------------------------------------------------------- */
     158                 : /*      Do we want to apply an approximating transformation?            */
     159                 : /* -------------------------------------------------------------------- */
     160              14 :     if( dfMaxError > 0.0 )
     161                 :     {
     162                 :         psWO->pTransformerArg = 
     163                 :             GDALCreateApproxTransformer( psWO->pfnTransformer, 
     164                 :                                          psWO->pTransformerArg, 
     165               2 :                                          dfMaxError );
     166               2 :         psWO->pfnTransformer = GDALApproxTransform;
     167               2 :         GDALApproxTransformerOwnsSubtransformer(psWO->pTransformerArg, TRUE);
     168                 :     }
     169                 : 
     170                 : /* -------------------------------------------------------------------- */
     171                 : /*      Create the VRT file.                                            */
     172                 : /* -------------------------------------------------------------------- */
     173                 :     GDALDatasetH hDstDS;
     174                 : 
     175                 :     hDstDS = GDALCreateWarpedVRT( hSrcDS, nDstPixels, nDstLines, 
     176              14 :                                   adfDstGeoTransform, psWO );
     177                 : 
     178              14 :     GDALDestroyWarpOptions( psWO );
     179                 : 
     180              14 :     if( pszDstWKT != NULL )
     181              13 :         GDALSetProjection( hDstDS, pszDstWKT );
     182               1 :     else if( pszSrcWKT != NULL )
     183               0 :         GDALSetProjection( hDstDS, pszDstWKT );
     184               1 :     else if( GDALGetGCPCount( hSrcDS ) > 0 )
     185               1 :         GDALSetProjection( hDstDS, GDALGetGCPProjection( hSrcDS ) );
     186                 :     else 
     187               0 :         GDALSetProjection( hDstDS, GDALGetProjectionRef( hSrcDS ) );
     188                 : 
     189              14 :     return hDstDS;
     190                 : }
     191                 : 
     192                 : /************************************************************************/
     193                 : /*                        GDALCreateWarpedVRT()                         */
     194                 : /************************************************************************/
     195                 : 
     196                 : /**
     197                 :  * Create virtual warped dataset.
     198                 :  *
     199                 :  * This function will create a warped virtual file representing the 
     200                 :  * input image warped based on a provided transformation.  Output bounds
     201                 :  * and resolution are provided explicitly.
     202                 :  *
     203                 :  * Note that the constructed GDALDatasetH will acquire one or more references 
     204                 :  * to the passed in hSrcDS.  Reference counting semantics on the source 
     205                 :  * dataset should be honoured.  That is, don't just GDALClose() it unless it 
     206                 :  * was opened with GDALOpenShared(). 
     207                 :  *
     208                 :  * @param hSrcDS The source dataset. 
     209                 :  *
     210                 :  * @param nPixels Width of the virtual warped dataset to create
     211                 :  *
     212                 :  * @param nLines Height of the virtual warped dataset to create
     213                 :  *
     214                 :  * @param padfGeoTransform Geotransform matrix of the virtual warped dataset to create
     215                 :  *
     216                 :  * @param psOptions Warp options. Must be different from NULL.
     217                 :  *
     218                 :  * @return NULL on failure, or a new virtual dataset handle on success.
     219                 :  */
     220                 : 
     221                 : GDALDatasetH CPL_STDCALL
     222              14 : GDALCreateWarpedVRT( GDALDatasetH hSrcDS, 
     223                 :                      int nPixels, int nLines, double *padfGeoTransform,
     224                 :                      GDALWarpOptions *psOptions )
     225                 :     
     226                 : {
     227              14 :     VALIDATE_POINTER1( hSrcDS, "GDALCreateWarpedVRT", NULL );
     228                 : 
     229                 : /* -------------------------------------------------------------------- */
     230                 : /*      Create the VRTDataset and populate it with bands.               */
     231                 : /* -------------------------------------------------------------------- */
     232              14 :     VRTWarpedDataset *poDS = new VRTWarpedDataset( nPixels, nLines );
     233                 :     int i;
     234                 : 
     235              14 :     psOptions->hDstDS = (GDALDatasetH) poDS;
     236                 : 
     237              14 :     poDS->SetGeoTransform( padfGeoTransform );
     238                 : 
     239              30 :     for( i = 0; i < psOptions->nBandCount; i++ )
     240                 :     {
     241                 :         VRTWarpedRasterBand *poBand;
     242                 :         GDALRasterBand *poSrcBand = (GDALRasterBand *) 
     243              16 :             GDALGetRasterBand( hSrcDS, i+1 );
     244                 : 
     245              16 :         poDS->AddBand( poSrcBand->GetRasterDataType(), NULL );
     246                 : 
     247              16 :         poBand = (VRTWarpedRasterBand *) poDS->GetRasterBand( i+1 );
     248              16 :         poBand->CopyCommonInfoFrom( poSrcBand );
     249                 :     }
     250                 : 
     251                 : /* -------------------------------------------------------------------- */
     252                 : /*      Initialize the warp on the VRTWarpedDataset.                    */
     253                 : /* -------------------------------------------------------------------- */
     254              14 :     poDS->Initialize( psOptions );
     255                 :     
     256              14 :     return (GDALDatasetH) poDS;
     257                 : }
     258                 : 
     259                 : /************************************************************************/
     260                 : /* ==================================================================== */
     261                 : /*                          VRTWarpedDataset                            */
     262                 : /* ==================================================================== */
     263                 : /************************************************************************/
     264                 : 
     265                 : /************************************************************************/
     266                 : /*                          VRTWarpedDataset()                          */
     267                 : /************************************************************************/
     268                 : 
     269              65 : VRTWarpedDataset::VRTWarpedDataset( int nXSize, int nYSize )
     270              65 :         : VRTDataset( nXSize, nYSize )
     271                 : 
     272                 : {
     273              65 :     poWarper = NULL;
     274              65 :     nBlockXSize = MIN(nXSize, 512);
     275              65 :     nBlockYSize = MIN(nYSize, 128);
     276              65 :     eAccess = GA_Update;
     277                 : 
     278              65 :     nOverviewCount = 0;
     279              65 :     papoOverviews = NULL;
     280              65 : }
     281                 : 
     282                 : /************************************************************************/
     283                 : /*                         ~VRTWarpedDataset()                          */
     284                 : /************************************************************************/
     285                 : 
     286              65 : VRTWarpedDataset::~VRTWarpedDataset()
     287                 : 
     288                 : {
     289              65 :     CloseDependentDatasets();
     290              65 : }
     291                 : 
     292                 : /************************************************************************/
     293                 : /*                        CloseDependentDatasets()                      */
     294                 : /************************************************************************/
     295                 : 
     296              65 : int VRTWarpedDataset::CloseDependentDatasets()
     297                 : {
     298              65 :     FlushCache();
     299                 : 
     300              65 :     int bHasDroppedRef = VRTDataset::CloseDependentDatasets();
     301                 : 
     302                 : /* -------------------------------------------------------------------- */
     303                 : /*      Cleanup overviews.                                              */
     304                 : /* -------------------------------------------------------------------- */
     305                 :     int iOverview;
     306                 : 
     307              67 :     for( iOverview = 0; iOverview < nOverviewCount; iOverview++ )
     308                 :     {
     309               2 :         GDALDatasetH hDS = (GDALDatasetH) papoOverviews[iOverview];
     310                 : 
     311               2 :         if( GDALDereferenceDataset( hDS ) < 1 )
     312                 :         {
     313               2 :             GDALReferenceDataset( hDS );
     314               2 :             GDALClose( hDS );
     315               2 :             bHasDroppedRef = TRUE;
     316                 :         }
     317                 :     }
     318                 : 
     319              65 :     CPLFree( papoOverviews );
     320              65 :     nOverviewCount = 0;
     321              65 :     papoOverviews = NULL;
     322                 : 
     323                 : /* -------------------------------------------------------------------- */
     324                 : /*      Cleanup warper if one is in effect.                             */
     325                 : /* -------------------------------------------------------------------- */
     326              65 :     if( poWarper != NULL )
     327                 :     {
     328              65 :         const GDALWarpOptions *psWO = poWarper->GetOptions();
     329                 : 
     330                 : /* -------------------------------------------------------------------- */
     331                 : /*      We take care to only call GDALClose() on psWO->hSrcDS if the    */
     332                 : /*      reference count drops to zero.  This is makes it so that we     */
     333                 : /*      can operate reference counting semantics more-or-less           */
     334                 : /*      properly even if the dataset isn't open in shared mode,         */
     335                 : /*      though we require that the caller also honour the reference     */
     336                 : /*      counting semantics even though it isn't a shared dataset.       */
     337                 : /* -------------------------------------------------------------------- */
     338              65 :         if( psWO->hSrcDS != NULL )
     339                 :         {
     340              65 :             if( GDALDereferenceDataset( psWO->hSrcDS ) < 1 )
     341                 :             {
     342              47 :                 GDALReferenceDataset( psWO->hSrcDS );
     343              47 :                 GDALClose( psWO->hSrcDS );
     344              47 :                 bHasDroppedRef = TRUE;
     345                 :             }
     346                 :         }
     347                 : 
     348                 : /* -------------------------------------------------------------------- */
     349                 : /*      We are responsible for cleaning up the transformer outselves.   */
     350                 : /* -------------------------------------------------------------------- */
     351              65 :         if( psWO->pTransformerArg != NULL )
     352              65 :             GDALDestroyTransformer( psWO->pTransformerArg );
     353                 : 
     354              65 :         delete poWarper;
     355              65 :         poWarper = NULL;
     356                 :     }
     357                 : 
     358                 : /* -------------------------------------------------------------------- */
     359                 : /*      Destroy the raster bands if they exist.                         */
     360                 : /* -------------------------------------------------------------------- */
     361             175 :     for( int iBand = 0; iBand < nBands; iBand++ )
     362                 :     {
     363             110 :        delete papoBands[iBand];
     364                 :     }
     365              65 :     nBands = 0;
     366                 : 
     367              65 :     return bHasDroppedRef;
     368                 : }
     369                 : 
     370                 : /************************************************************************/
     371                 : /*                             Initialize()                             */
     372                 : /*                                                                      */
     373                 : /*      Initialize a dataset from passed in warp options.               */
     374                 : /************************************************************************/
     375                 : 
     376              20 : CPLErr VRTWarpedDataset::Initialize( void *psWO )
     377                 : 
     378                 : {
     379              20 :     if( poWarper != NULL )
     380               0 :         delete poWarper;
     381                 : 
     382              20 :     poWarper = new GDALWarpOperation();
     383                 : 
     384              20 :     GDALWarpOptions* psWO_Dup = GDALCloneWarpOptions((GDALWarpOptions *) psWO);
     385                 : 
     386                 :     /* Avoid errors when adding an alpha band, but source dataset has */
     387                 :     /* no alpha band (#4571) */
     388              20 :     if (CSLFetchNameValue( psWO_Dup->papszWarpOptions, "INIT_DEST" ) == NULL)
     389              14 :         psWO_Dup->papszWarpOptions = CSLSetNameValue(psWO_Dup->papszWarpOptions, "INIT_DEST", "0");
     390                 : 
     391                 :     // The act of initializing this warped dataset with this warp options
     392                 :     // will result in our assuming ownership of a reference to the
     393                 :     // hSrcDS.
     394                 : 
     395              20 :     if( ((GDALWarpOptions *) psWO)->hSrcDS != NULL )
     396              20 :         GDALReferenceDataset( psWO_Dup->hSrcDS );
     397                 : 
     398              20 :     CPLErr eErr = poWarper->Initialize( psWO_Dup );
     399                 : 
     400              20 :     GDALDestroyWarpOptions(psWO_Dup);
     401                 : 
     402              20 :     return eErr;
     403                 : }
     404                 : 
     405                 : /************************************************************************/
     406                 : /*                            GetFileList()                             */
     407                 : /************************************************************************/
     408                 : 
     409               0 : char** VRTWarpedDataset::GetFileList()
     410                 : {
     411               0 :     char** papszFileList = GDALDataset::GetFileList();
     412                 :     
     413               0 :     if( poWarper != NULL )
     414                 :     {
     415               0 :         const GDALWarpOptions *psWO = poWarper->GetOptions();
     416                 :         const char* pszFilename;
     417                 :         
     418               0 :         if( psWO->hSrcDS != NULL &&
     419                 :             (pszFilename =
     420               0 :                     ((GDALDataset*)psWO->hSrcDS)->GetDescription()) != NULL )
     421                 :         {
     422                 :             VSIStatBufL  sStat;
     423               0 :             if( VSIStatL( pszFilename, &sStat ) == 0 )
     424                 :             {
     425               0 :                 papszFileList = CSLAddString(papszFileList, pszFilename);
     426                 :             }
     427                 :         }
     428                 :     }
     429                 :     
     430               0 :     return papszFileList;
     431                 : }
     432                 : 
     433                 : 
     434                 : 
     435                 : /************************************************************************/
     436                 : /* ==================================================================== */
     437                 : /*                    VRTWarpedOverviewTransformer                      */
     438                 : /* ==================================================================== */
     439                 : /************************************************************************/
     440                 : 
     441                 : typedef struct {
     442                 :     GDALTransformerInfo sTI;
     443                 : 
     444                 :     GDALTransformerFunc pfnBaseTransformer;
     445                 :     void              *pBaseTransformerArg;
     446                 :     int                bOwnSubtransformer;
     447                 : 
     448                 :     double            dfXOverviewFactor;
     449                 :     double            dfYOverviewFactor;
     450                 : } VWOTInfo;
     451                 : 
     452                 : 
     453                 : static
     454                 : void* VRTCreateWarpedOverviewTransformer( GDALTransformerFunc pfnBaseTransformer,
     455                 :                                           void *pBaseTransformArg,
     456                 :                                           double dfXOverviewFactor,
     457                 :                                           double dfYOverviewFactor );
     458                 : static
     459                 : void VRTDestroyWarpedOverviewTransformer(void* pTransformArg);
     460                 : 
     461                 : /************************************************************************/
     462                 : /*                VRTSerializeWarpedOverviewTransformer()               */
     463                 : /************************************************************************/
     464                 : 
     465                 : static CPLXMLNode *
     466               0 : VRTSerializeWarpedOverviewTransformer( void *pTransformArg )
     467                 : 
     468                 : {
     469                 :     CPLXMLNode *psTree;
     470               0 :     VWOTInfo *psInfo = (VWOTInfo *) pTransformArg;
     471                 : 
     472               0 :     psTree = CPLCreateXMLNode( NULL, CXT_Element, "WarpedOverviewTransformer" );
     473                 : 
     474                 :     CPLCreateXMLElementAndValue( psTree, "XFactor",
     475               0 :                                  CPLString().Printf("%g",psInfo->dfXOverviewFactor) );
     476                 :     CPLCreateXMLElementAndValue( psTree, "YFactor",
     477               0 :                                  CPLString().Printf("%g",psInfo->dfYOverviewFactor) );
     478                 : 
     479                 : /* -------------------------------------------------------------------- */
     480                 : /*      Capture underlying transformer.                                 */
     481                 : /* -------------------------------------------------------------------- */
     482                 :     CPLXMLNode *psTransformerContainer;
     483                 :     CPLXMLNode *psTransformer;
     484                 : 
     485                 :     psTransformerContainer =
     486               0 :         CPLCreateXMLNode( psTree, CXT_Element, "BaseTransformer" );
     487                 : 
     488                 :     psTransformer = GDALSerializeTransformer( psInfo->pfnBaseTransformer,
     489               0 :                                               psInfo->pBaseTransformerArg );
     490               0 :     if( psTransformer != NULL )
     491               0 :         CPLAddXMLChild( psTransformerContainer, psTransformer );
     492                 : 
     493               0 :     return psTree;
     494                 : }
     495                 : 
     496                 : /************************************************************************/
     497                 : /*           VRTWarpedOverviewTransformerOwnsSubtransformer()           */
     498                 : /************************************************************************/
     499                 : 
     500               0 : static void VRTWarpedOverviewTransformerOwnsSubtransformer( void *pTransformArg,
     501                 :                                                           int bOwnFlag )
     502                 : {
     503                 :     VWOTInfo *psInfo =
     504               0 :             (VWOTInfo *) pTransformArg;
     505                 : 
     506               0 :     psInfo->bOwnSubtransformer = bOwnFlag;
     507               0 : }
     508                 : 
     509                 : /************************************************************************/
     510                 : /*            VRTDeserializeWarpedOverviewTransformer()                 */
     511                 : /************************************************************************/
     512                 : 
     513               0 : void* VRTDeserializeWarpedOverviewTransformer( CPLXMLNode *psTree )
     514                 : 
     515                 : {
     516               0 :     double dfXOverviewFactor = atof(CPLGetXMLValue( psTree, "XFactor",  "1" ));
     517               0 :     double dfYOverviewFactor = atof(CPLGetXMLValue( psTree, "YFactor",  "1" ));
     518                 :     CPLXMLNode *psContainer;
     519               0 :     GDALTransformerFunc pfnBaseTransform = NULL;
     520               0 :     void *pBaseTransformerArg = NULL;
     521                 : 
     522               0 :     psContainer = CPLGetXMLNode( psTree, "BaseTransformer" );
     523                 : 
     524               0 :     if( psContainer != NULL && psContainer->psChild != NULL )
     525                 :     {
     526                 :         GDALDeserializeTransformer( psContainer->psChild,
     527                 :                                     &pfnBaseTransform,
     528               0 :                                     &pBaseTransformerArg );
     529                 : 
     530                 :     }
     531                 : 
     532               0 :     if( pfnBaseTransform == NULL )
     533                 :     {
     534                 :         CPLError( CE_Failure, CPLE_AppDefined,
     535               0 :                   "Cannot get base transform for scaled coord transformer." );
     536               0 :         return NULL;
     537                 :     }
     538                 :     else
     539                 :     {
     540                 :         void *pApproxCBData =
     541                 :                        VRTCreateWarpedOverviewTransformer( pfnBaseTransform,
     542                 :                                                            pBaseTransformerArg,
     543                 :                                                            dfXOverviewFactor,
     544               0 :                                                            dfYOverviewFactor );
     545               0 :         VRTWarpedOverviewTransformerOwnsSubtransformer( pApproxCBData, TRUE );
     546                 : 
     547               0 :         return pApproxCBData;
     548                 :     }
     549                 : }
     550                 : 
     551                 : 
     552                 : /************************************************************************/
     553                 : /*                   VRTCreateWarpedOverviewTransformer()               */
     554                 : /************************************************************************/
     555                 : 
     556                 : static
     557               2 : void* VRTCreateWarpedOverviewTransformer( GDALTransformerFunc pfnBaseTransformer,
     558                 :                                           void *pBaseTransformerArg,
     559                 :                                           double dfXOverviewFactor,
     560                 :                                           double dfYOverviewFactor)
     561                 : 
     562                 : {
     563                 :     VWOTInfo *psSCTInfo;
     564                 : 
     565               2 :     if (pfnBaseTransformer == NULL)
     566               0 :         return NULL;
     567                 : 
     568                 :     psSCTInfo = (VWOTInfo*)
     569               2 :                     CPLMalloc(sizeof(VWOTInfo));
     570               2 :     psSCTInfo->pfnBaseTransformer = pfnBaseTransformer;
     571               2 :     psSCTInfo->pBaseTransformerArg = pBaseTransformerArg;
     572               2 :     psSCTInfo->dfXOverviewFactor = dfXOverviewFactor;
     573               2 :     psSCTInfo->dfYOverviewFactor = dfYOverviewFactor;
     574               2 :     psSCTInfo->bOwnSubtransformer = FALSE;
     575                 : 
     576               2 :     strcpy( psSCTInfo->sTI.szSignature, "GTI" );
     577               2 :     psSCTInfo->sTI.pszClassName = "VRTWarpedOverviewTransformer";
     578               2 :     psSCTInfo->sTI.pfnTransform = VRTWarpedOverviewTransform;
     579               2 :     psSCTInfo->sTI.pfnCleanup = VRTDestroyWarpedOverviewTransformer;
     580               2 :     psSCTInfo->sTI.pfnSerialize = VRTSerializeWarpedOverviewTransformer;
     581                 : 
     582               2 :     return psSCTInfo;
     583                 : }
     584                 : 
     585                 : /************************************************************************/
     586                 : /*               VRTDestroyWarpedOverviewTransformer()                  */
     587                 : /************************************************************************/
     588                 : 
     589                 : static
     590               2 : void VRTDestroyWarpedOverviewTransformer(void* pTransformArg)
     591                 : {
     592               2 :     VWOTInfo *psInfo = (VWOTInfo *) pTransformArg;
     593                 : 
     594               2 :     if( psInfo->bOwnSubtransformer )
     595               0 :         GDALDestroyTransformer( psInfo->pBaseTransformerArg );
     596                 : 
     597               2 :     CPLFree( psInfo );
     598               2 : }
     599                 : 
     600                 : /************************************************************************/
     601                 : /*                     VRTWarpedOverviewTransform()                     */
     602                 : /************************************************************************/
     603                 : 
     604             258 : int VRTWarpedOverviewTransform( void *pTransformArg, int bDstToSrc,
     605                 :                                 int nPointCount,
     606                 :                                 double *padfX, double *padfY, double *padfZ,
     607                 :                                 int *panSuccess )
     608                 : 
     609                 : {
     610             258 :     VWOTInfo *psInfo = (VWOTInfo *) pTransformArg;
     611                 :     int i, bSuccess;
     612                 : 
     613             258 :     if( bDstToSrc )
     614                 :     {
     615           64426 :         for( i = 0; i < nPointCount; i++ )
     616                 :         {
     617           64168 :             padfX[i] *= psInfo->dfXOverviewFactor;
     618           64168 :             padfY[i] *= psInfo->dfYOverviewFactor;
     619                 :         }
     620                 :     }
     621                 : 
     622                 :     bSuccess = psInfo->pfnBaseTransformer( psInfo->pBaseTransformerArg,
     623                 :                                            bDstToSrc,
     624                 :                                            nPointCount, padfX, padfY, padfZ,
     625             258 :                                            panSuccess );
     626                 : 
     627             258 :     if( !bDstToSrc )
     628                 :     {
     629               0 :         for( i = 0; i < nPointCount; i++ )
     630                 :         {
     631               0 :             padfX[i] /= psInfo->dfXOverviewFactor;
     632               0 :             padfY[i] /= psInfo->dfYOverviewFactor;
     633                 :         }
     634                 :     }
     635                 : 
     636             258 :     return bSuccess;
     637                 : }
     638                 : 
     639                 : /************************************************************************/
     640                 : /*                           BuildOverviews()                           */
     641                 : /*                                                                      */
     642                 : /*      For overviews, we actually just build a whole new dataset       */
     643                 : /*      with an extra layer of transformation on the warper used to     */
     644                 : /*      accomplish downsampling by the desired factor.                  */
     645                 : /************************************************************************/
     646                 : 
     647                 : CPLErr 
     648               2 : VRTWarpedDataset::IBuildOverviews( const char *pszResampling, 
     649                 :                                    int nOverviews, int *panOverviewList, 
     650                 :                                    int nListBands, int *panBandList,
     651                 :                                    GDALProgressFunc pfnProgress, 
     652                 :                                    void * pProgressData )
     653                 :     
     654                 : {
     655                 : /* -------------------------------------------------------------------- */
     656                 : /*      Initial progress result.                                        */
     657                 : /* -------------------------------------------------------------------- */
     658               2 :     if( !pfnProgress( 0.0, NULL, pProgressData ) )
     659                 :     {
     660               0 :         CPLError( CE_Failure, CPLE_UserInterrupt, "User terminated" );
     661               0 :         return CE_Failure;
     662                 :     }
     663                 : 
     664                 : /* -------------------------------------------------------------------- */
     665                 : /*      Establish which of the overview levels we already have, and     */
     666                 : /*      which are new.                                                  */
     667                 : /* -------------------------------------------------------------------- */
     668               2 :     int   i, nNewOverviews, *panNewOverviewList = NULL;
     669                 : 
     670               2 :     nNewOverviews = 0;
     671               2 :     panNewOverviewList = (int *) CPLCalloc(sizeof(int),nOverviews);
     672               4 :     for( i = 0; i < nOverviews; i++ )
     673                 :     {
     674                 :         int   j;
     675                 : 
     676               2 :         for( j = 0; j < nOverviewCount; j++ )
     677                 :         {
     678                 :             int    nOvFactor;
     679               0 :             VRTWarpedDataset *poOverview = papoOverviews[j];
     680                 :             
     681                 :             nOvFactor = (int) 
     682               0 :                 (0.5+GetRasterXSize() / (double) poOverview->GetRasterXSize());
     683                 : 
     684               0 :             if( nOvFactor == panOverviewList[i] 
     685                 :                 || nOvFactor == GDALOvLevelAdjust( panOverviewList[i], 
     686                 :                                                    GetRasterXSize() ) )
     687               0 :                 panOverviewList[i] *= -1;
     688                 :         }
     689                 : 
     690               2 :         if( panOverviewList[i] > 0 )
     691               2 :             panNewOverviewList[nNewOverviews++] = panOverviewList[i];
     692                 :     }
     693                 : 
     694                 : /* -------------------------------------------------------------------- */
     695                 : /*      Create each missing overview (we don't need to do anything      */
     696                 : /*      to update existing overviews).                                  */
     697                 : /* -------------------------------------------------------------------- */
     698               4 :     for( i = 0; i < nNewOverviews; i++ )
     699                 :     {
     700                 :         int    nOXSize, nOYSize, iBand;
     701                 :         VRTWarpedDataset *poOverviewDS;
     702                 :         
     703                 : /* -------------------------------------------------------------------- */
     704                 : /*      What size should this overview be.                              */
     705                 : /* -------------------------------------------------------------------- */
     706               2 :         nOXSize = (GetRasterXSize() + panNewOverviewList[i] - 1) 
     707               4 :             / panNewOverviewList[i];
     708                 :                                  
     709               2 :         nOYSize = (GetRasterYSize() + panNewOverviewList[i] - 1) 
     710               4 :             / panNewOverviewList[i];
     711                 : 
     712                 : /* -------------------------------------------------------------------- */
     713                 : /*      Create the overview dataset.                                    */
     714                 : /* -------------------------------------------------------------------- */
     715               2 :         poOverviewDS = new VRTWarpedDataset( nOXSize, nOYSize );
     716                 :         
     717               8 :         for( iBand = 0; iBand < GetRasterCount(); iBand++ )
     718                 :         {
     719               2 :             GDALRasterBand *poOldBand = GetRasterBand(iBand+1);
     720                 :             VRTWarpedRasterBand *poNewBand = 
     721                 :                 new VRTWarpedRasterBand( poOverviewDS, iBand+1, 
     722               2 :                                          poOldBand->GetRasterDataType() );
     723                 :             
     724               2 :             poNewBand->CopyCommonInfoFrom( poOldBand );
     725               2 :             poOverviewDS->SetBand( iBand+1, poNewBand );
     726                 :         }
     727                 : 
     728               2 :         nOverviewCount++;
     729                 :         papoOverviews = (VRTWarpedDataset **)
     730               2 :             CPLRealloc( papoOverviews, sizeof(void*) * nOverviewCount );
     731                 : 
     732               2 :         papoOverviews[nOverviewCount-1] = poOverviewDS;
     733                 :         
     734                 : /* -------------------------------------------------------------------- */
     735                 : /*      Prepare update transformation information that will apply       */
     736                 : /*      the overview decimation.                                        */
     737                 : /* -------------------------------------------------------------------- */
     738               2 :         GDALWarpOptions *psWO = (GDALWarpOptions *) poWarper->GetOptions();
     739                 : 
     740                 : /* -------------------------------------------------------------------- */
     741                 : /*      Initialize the new dataset with adjusted warp options, and      */
     742                 : /*      then restore to original condition.                             */
     743                 : /* -------------------------------------------------------------------- */
     744               2 :         GDALTransformerFunc pfnTransformerBase = psWO->pfnTransformer;
     745               2 :         void* pTransformerBaseArg = psWO->pTransformerArg;
     746                 : 
     747               2 :         psWO->pfnTransformer = VRTWarpedOverviewTransform;
     748                 :         psWO->pTransformerArg = VRTCreateWarpedOverviewTransformer(
     749                 :                                         pfnTransformerBase,
     750                 :                                         pTransformerBaseArg,
     751                 :                                         GetRasterXSize() / (double) nOXSize,
     752               2 :                                         GetRasterYSize() / (double) nOYSize );
     753                 : 
     754               2 :         poOverviewDS->Initialize( psWO );
     755                 : 
     756               2 :         psWO->pfnTransformer = pfnTransformerBase;
     757               2 :         psWO->pTransformerArg = pTransformerBaseArg;
     758                 :     }
     759                 : 
     760               2 :     CPLFree( panNewOverviewList );
     761                 : 
     762                 : /* -------------------------------------------------------------------- */
     763                 : /*      Progress finished.                                              */
     764                 : /* -------------------------------------------------------------------- */
     765               2 :     pfnProgress( 1.0, NULL, pProgressData );
     766                 : 
     767               2 :     SetNeedsFlush();
     768                 : 
     769               2 :     return CE_None;
     770                 : }
     771                 : 
     772                 : /************************************************************************/
     773                 : /*                      GDALInitializeWarpedVRT()                       */
     774                 : /************************************************************************/
     775                 : 
     776                 : /**
     777                 :  * Set warp info on virtual warped dataset.
     778                 :  *
     779                 :  * Initializes all the warping information for a virtual warped dataset.
     780                 :  *
     781                 :  * This method is the same as the C++ method VRTWarpedDataset::Initialize().
     782                 :  *
     783                 :  * @param hDS dataset previously created with the VRT driver, and a 
     784                 :  * SUBCLASS of "VRTWarpedDataset".
     785                 :  * 
     786                 :  * @param psWO the warp options to apply.  Note that ownership of the
     787                 :  * transformation information is taken over by the function though everything
     788                 :  * else remains the property of the caller. 
     789                 :  *
     790                 :  * @return CE_None on success or CE_Failure if an error occurs. 
     791                 :  */
     792                 : 
     793                 : CPLErr CPL_STDCALL 
     794               4 : GDALInitializeWarpedVRT( GDALDatasetH hDS, GDALWarpOptions *psWO )
     795                 : 
     796                 : {
     797               4 :     VALIDATE_POINTER1( hDS, "GDALInitializeWarpedVRT", CE_Failure );
     798                 : 
     799               4 :     return ((VRTWarpedDataset *) hDS)->Initialize( psWO );
     800                 : }
     801                 : 
     802                 : /************************************************************************/
     803                 : /*                              XMLInit()                               */
     804                 : /************************************************************************/
     805                 : 
     806              45 : CPLErr VRTWarpedDataset::XMLInit( CPLXMLNode *psTree, const char *pszVRTPath )
     807                 : 
     808                 : {
     809                 :     CPLErr eErr;
     810                 : 
     811                 : /* -------------------------------------------------------------------- */
     812                 : /*      Initialize blocksize before calling sub-init so that the        */
     813                 : /*      band initializers can get it from the dataset object when       */
     814                 : /*      they are created.                                               */
     815                 : /* -------------------------------------------------------------------- */
     816              45 :     nBlockXSize = atoi(CPLGetXMLValue(psTree,"BlockXSize","512"));
     817              45 :     nBlockYSize = atoi(CPLGetXMLValue(psTree,"BlockYSize","128"));
     818                 : 
     819                 : /* -------------------------------------------------------------------- */
     820                 : /*      Initialize all the general VRT stuff.  This will even           */
     821                 : /*      create the VRTWarpedRasterBands and initialize them.            */
     822                 : /* -------------------------------------------------------------------- */
     823              45 :     eErr = VRTDataset::XMLInit( psTree, pszVRTPath );
     824                 : 
     825              45 :     if( eErr != CE_None )
     826               0 :         return eErr;
     827                 : 
     828                 : /* -------------------------------------------------------------------- */
     829                 : /*      Find the GDALWarpOptions XML tree.                              */
     830                 : /* -------------------------------------------------------------------- */
     831                 :     CPLXMLNode *psOptionsTree;
     832              45 :     psOptionsTree = CPLGetXMLNode( psTree, "GDALWarpOptions" );
     833              45 :     if( psOptionsTree == NULL )
     834                 :     {
     835                 :         CPLError( CE_Failure, CPLE_AppDefined,
     836               0 :                   "Count not find required GDALWarpOptions in XML." );
     837               0 :         return CE_Failure;
     838                 :     }
     839                 : 
     840                 : /* -------------------------------------------------------------------- */
     841                 : /*      Adjust the SourceDataset in the warp options to take into       */
     842                 : /*      account that it is relative to the VRT if appropriate.          */
     843                 : /* -------------------------------------------------------------------- */
     844                 :     int bRelativeToVRT = 
     845                 :         atoi(CPLGetXMLValue(psOptionsTree,
     846              45 :                             "SourceDataset.relativeToVRT", "0" ));
     847                 : 
     848                 :     const char *pszRelativePath = CPLGetXMLValue(psOptionsTree,
     849              45 :                                                  "SourceDataset", "" );
     850                 :     char *pszAbsolutePath;
     851                 : 
     852              45 :     if( bRelativeToVRT )
     853                 :         pszAbsolutePath = 
     854                 :             CPLStrdup(CPLProjectRelativeFilename( pszVRTPath, 
     855              42 :                                                   pszRelativePath ) );
     856                 :     else
     857               3 :         pszAbsolutePath = CPLStrdup(pszRelativePath);
     858                 : 
     859              45 :     CPLSetXMLValue( psOptionsTree, "SourceDataset", pszAbsolutePath );
     860              45 :     CPLFree( pszAbsolutePath );
     861                 : 
     862                 : /* -------------------------------------------------------------------- */
     863                 : /*      And instantiate the warp options, and corresponding warp        */
     864                 : /*      operation.                                                      */
     865                 : /* -------------------------------------------------------------------- */
     866                 :     GDALWarpOptions *psWO;
     867                 : 
     868              45 :     psWO = GDALDeserializeWarpOptions( psOptionsTree );
     869              45 :     if( psWO == NULL )
     870               0 :         return CE_Failure;
     871                 : 
     872                 :     /* Avoid errors when adding an alpha band, but source dataset has */
     873                 :     /* no alpha band (#4571) */
     874              45 :     if (CSLFetchNameValue( psWO->papszWarpOptions, "INIT_DEST" ) == NULL)
     875              18 :         psWO->papszWarpOptions = CSLSetNameValue(psWO->papszWarpOptions, "INIT_DEST", "0");
     876                 : 
     877              45 :     this->eAccess = GA_Update;
     878                 : 
     879              45 :     if( psWO->hDstDS != NULL )
     880                 :     {
     881               0 :         GDALClose( psWO->hDstDS );
     882               0 :         psWO->hDstDS = NULL;
     883                 :     }
     884                 : 
     885              45 :     psWO->hDstDS = this;
     886                 : 
     887                 : /* -------------------------------------------------------------------- */
     888                 : /*      Instantiate the warp operation.                                 */
     889                 : /* -------------------------------------------------------------------- */
     890              45 :     poWarper = new GDALWarpOperation();
     891                 : 
     892              45 :     eErr = poWarper->Initialize( psWO );
     893              45 :     if( eErr != CE_None)
     894                 :     {
     895                 : /* -------------------------------------------------------------------- */
     896                 : /*      We are responsible for cleaning up the transformer outselves.   */
     897                 : /* -------------------------------------------------------------------- */
     898               0 :         if( psWO->pTransformerArg != NULL )
     899                 :         {
     900               0 :             GDALDestroyTransformer( psWO->pTransformerArg );
     901               0 :             psWO->pTransformerArg = NULL;
     902                 :         }
     903                 : 
     904               0 :         if( psWO->hSrcDS != NULL )
     905                 :         {
     906               0 :             GDALClose( psWO->hSrcDS );
     907               0 :             psWO->hSrcDS = NULL;
     908                 :         }
     909                 :     }
     910                 : 
     911              45 :     GDALDestroyWarpOptions( psWO );
     912              45 :     if( eErr != CE_None )
     913                 :     {
     914               0 :         delete poWarper;
     915               0 :         poWarper = NULL;
     916                 :     }
     917                 : 
     918                 : /* -------------------------------------------------------------------- */
     919                 : /*      Generate overviews, if appropriate.                             */
     920                 : /* -------------------------------------------------------------------- */
     921                 :     char **papszTokens = CSLTokenizeString( 
     922              45 :         CPLGetXMLValue( psTree, "OverviewList", "" ));
     923                 :     int iOverview;
     924                 : 
     925              92 :     for( iOverview = 0; 
     926              46 :          papszTokens != NULL && papszTokens[iOverview] != NULL;
     927                 :          iOverview++ )
     928                 :     {
     929               1 :         int nOvFactor = atoi(papszTokens[iOverview]);
     930                 : 
     931               1 :         if (nOvFactor > 0)
     932               1 :             BuildOverviews( "NEAREST", 1, &nOvFactor, 0, NULL, NULL, NULL );
     933                 :         else
     934                 :             CPLError(CE_Failure, CPLE_AppDefined,
     935               0 :                      "Bad value for overview factor : %s", papszTokens[iOverview]);
     936                 :     }
     937                 : 
     938              45 :     CSLDestroy( papszTokens );
     939                 : 
     940              45 :     return eErr;
     941                 : }
     942                 : 
     943                 : /************************************************************************/
     944                 : /*                           SerializeToXML()                           */
     945                 : /************************************************************************/
     946                 : 
     947               6 : CPLXMLNode *VRTWarpedDataset::SerializeToXML( const char *pszVRTPath )
     948                 : 
     949                 : {
     950                 :     CPLXMLNode *psTree;
     951                 : 
     952               6 :     psTree = VRTDataset::SerializeToXML( pszVRTPath );
     953                 : 
     954               6 :     if( psTree == NULL )
     955               0 :         return psTree;
     956                 : 
     957                 : /* -------------------------------------------------------------------- */
     958                 : /*      Set subclass.                                                   */
     959                 : /* -------------------------------------------------------------------- */
     960                 :     CPLCreateXMLNode( 
     961                 :         CPLCreateXMLNode( psTree, CXT_Attribute, "subClass" ), 
     962               6 :         CXT_Text, "VRTWarpedDataset" );
     963                 : 
     964                 : /* -------------------------------------------------------------------- */
     965                 : /*      Serialize the block size.                                       */
     966                 : /* -------------------------------------------------------------------- */
     967                 :     CPLCreateXMLElementAndValue( psTree, "BlockXSize",
     968               6 :                                  CPLSPrintf( "%d", nBlockXSize ) );
     969                 :     CPLCreateXMLElementAndValue( psTree, "BlockYSize",
     970               6 :                                  CPLSPrintf( "%d", nBlockYSize ) );
     971                 : 
     972                 : /* -------------------------------------------------------------------- */
     973                 : /*      Serialize the overview list.                                    */
     974                 : /* -------------------------------------------------------------------- */
     975               6 :     if( nOverviewCount > 0 )
     976                 :     {
     977                 :         char *pszOverviewList;
     978                 :         int iOverview;
     979                 :         
     980               1 :         pszOverviewList = (char *) CPLMalloc(nOverviewCount*8 + 10);
     981               1 :         pszOverviewList[0] = '\0';
     982               2 :         for( iOverview = 0; iOverview < nOverviewCount; iOverview++ )
     983                 :         {
     984                 :             int nOvFactor;
     985                 : 
     986                 :             nOvFactor = (int) 
     987                 :                 (0.5+GetRasterXSize() 
     988               1 :                  / (double) papoOverviews[iOverview]->GetRasterXSize());
     989                 : 
     990                 :             sprintf( pszOverviewList + strlen(pszOverviewList), 
     991               1 :                      "%d ", nOvFactor );
     992                 :         }
     993                 : 
     994               1 :         CPLCreateXMLElementAndValue( psTree, "OverviewList", pszOverviewList );
     995                 : 
     996               1 :         CPLFree( pszOverviewList );
     997                 :     }
     998                 : 
     999                 : /* ==================================================================== */
    1000                 : /*      Serialize the warp options.                                     */
    1001                 : /* ==================================================================== */
    1002                 :     CPLXMLNode *psWOTree;
    1003                 : 
    1004               6 :     if( poWarper != NULL )
    1005                 :     {
    1006                 : /* -------------------------------------------------------------------- */
    1007                 : /*      We reset the destination dataset name so it doesn't get         */
    1008                 : /*      written out in the serialize warp options.                      */
    1009                 : /* -------------------------------------------------------------------- */
    1010               6 :         char *pszSavedName = CPLStrdup(GetDescription());
    1011               6 :         SetDescription("");
    1012                 : 
    1013               6 :         psWOTree = GDALSerializeWarpOptions( poWarper->GetOptions() );
    1014               6 :         CPLAddXMLChild( psTree, psWOTree );
    1015                 : 
    1016               6 :         SetDescription( pszSavedName );
    1017               6 :         CPLFree( pszSavedName );
    1018                 : 
    1019                 : /* -------------------------------------------------------------------- */
    1020                 : /*      We need to consider making the source dataset relative to       */
    1021                 : /*      the VRT file if possible.  Adjust accordingly.                  */
    1022                 : /* -------------------------------------------------------------------- */
    1023               6 :         CPLXMLNode *psSDS = CPLGetXMLNode( psWOTree, "SourceDataset" );
    1024                 :         int bRelativeToVRT;
    1025                 :         char *pszRelativePath;
    1026                 : 
    1027                 :         pszRelativePath = 
    1028                 :             CPLStrdup(
    1029                 :                 CPLExtractRelativePath( pszVRTPath, psSDS->psChild->pszValue, 
    1030               6 :                                         &bRelativeToVRT ) );
    1031                 : 
    1032               6 :         CPLFree( psSDS->psChild->pszValue );
    1033               6 :         psSDS->psChild->pszValue = pszRelativePath;
    1034                 : 
    1035                 :         CPLCreateXMLNode( 
    1036                 :             CPLCreateXMLNode( psSDS, CXT_Attribute, "relativeToVRT" ), 
    1037               6 :             CXT_Text, bRelativeToVRT ? "1" : "0" );
    1038                 :     }
    1039                 : 
    1040               6 :     return psTree;
    1041                 : }
    1042                 : 
    1043                 : /************************************************************************/
    1044                 : /*                            GetBlockSize()                            */
    1045                 : /************************************************************************/
    1046                 : 
    1047             110 : void VRTWarpedDataset::GetBlockSize( int *pnBlockXSize, int *pnBlockYSize )
    1048                 : 
    1049                 : {
    1050             110 :     assert( NULL != pnBlockXSize );
    1051             110 :     assert( NULL != pnBlockYSize );
    1052                 : 
    1053             110 :     *pnBlockXSize = nBlockXSize;
    1054             110 :     *pnBlockYSize = nBlockYSize;
    1055             110 : }
    1056                 : 
    1057                 : /************************************************************************/
    1058                 : /*                            ProcessBlock()                            */
    1059                 : /*                                                                      */
    1060                 : /*      Warp a single requested block, and then push each band of       */
    1061                 : /*      the result into the block cache.                                */
    1062                 : /************************************************************************/
    1063                 : 
    1064              99 : CPLErr VRTWarpedDataset::ProcessBlock( int iBlockX, int iBlockY )
    1065                 : 
    1066                 : {
    1067              99 :     if( poWarper == NULL )
    1068               0 :         return CE_Failure;
    1069                 : 
    1070              99 :     const GDALWarpOptions *psWO = poWarper->GetOptions();
    1071                 : 
    1072                 : /* -------------------------------------------------------------------- */
    1073                 : /*      Allocate block of memory large enough to hold all the bands     */
    1074                 : /*      for this block.                                                 */
    1075                 : /* -------------------------------------------------------------------- */
    1076                 :     int iBand;
    1077                 :     GByte *pabyDstBuffer;
    1078                 :     int   nDstBufferSize;
    1079              99 :     int   nWordSize = (GDALGetDataTypeSize(psWO->eWorkingDataType) / 8);
    1080                 : 
    1081                 :     // FIXME? : risk of overflow in multiplication if nBlockXSize or nBlockYSize are very large
    1082              99 :     nDstBufferSize = nBlockXSize * nBlockYSize * psWO->nBandCount * nWordSize;
    1083                 : 
    1084              99 :     pabyDstBuffer = (GByte *) VSIMalloc(nDstBufferSize);
    1085                 : 
    1086              99 :     if( pabyDstBuffer == NULL )
    1087                 :     {
    1088                 :         CPLError( CE_Failure, CPLE_OutOfMemory,
    1089                 :                   "Out of memory allocating %d byte buffer in VRTWarpedDataset::ProcessBlock()",
    1090               0 :                   nDstBufferSize );
    1091               0 :         return CE_Failure;
    1092                 :     }       
    1093                 : 
    1094              99 :     memset( pabyDstBuffer, 0, nDstBufferSize );
    1095                 : 
    1096                 : /* -------------------------------------------------------------------- */
    1097                 : /*      Process INIT_DEST option to initialize the buffer prior to      */
    1098                 : /*      warping into it.                                                */
    1099                 : /* NOTE:The following code is 99% similar in gdalwarpoperation.cpp and  */
    1100                 : /*      vrtwarped.cpp. Be careful to keep it in sync !                  */
    1101                 : /* -------------------------------------------------------------------- */
    1102                 :     const char *pszInitDest = CSLFetchNameValue( psWO->papszWarpOptions,
    1103              99 :                                                  "INIT_DEST" );
    1104                 : 
    1105              99 :     if( pszInitDest != NULL && !EQUAL(pszInitDest, "") )
    1106                 :     {
    1107                 :         char **papszInitValues = 
    1108              99 :             CSLTokenizeStringComplex( pszInitDest, ",", FALSE, FALSE );
    1109              99 :         int nInitCount = CSLCount(papszInitValues);
    1110                 :                                                            
    1111             227 :         for( iBand = 0; iBand < psWO->nBandCount; iBand++ )
    1112                 :         {
    1113                 :             double adfInitRealImag[2];
    1114                 :             GByte *pBandData;
    1115             128 :             int nBandSize = nBlockXSize * nBlockYSize * nWordSize;
    1116             128 :             const char *pszBandInit = papszInitValues[MIN(iBand,nInitCount-1)];
    1117                 : 
    1118             136 :             if( EQUAL(pszBandInit,"NO_DATA")
    1119                 :                 && psWO->padfDstNoDataReal != NULL )
    1120                 :             {
    1121               8 :                 adfInitRealImag[0] = psWO->padfDstNoDataReal[iBand];
    1122               8 :                 adfInitRealImag[1] = psWO->padfDstNoDataImag[iBand];
    1123                 :             }
    1124                 :             else
    1125                 :             {
    1126                 :                 CPLStringToComplex( pszBandInit,
    1127             120 :                                     adfInitRealImag + 0, adfInitRealImag + 1);
    1128                 :             }
    1129                 : 
    1130             128 :             pBandData = ((GByte *) pabyDstBuffer) + iBand * nBandSize;
    1131                 :             
    1132             128 :             if( psWO->eWorkingDataType == GDT_Byte )
    1133                 :                 memset( pBandData, 
    1134             170 :                         MAX(0,MIN(255,(int)adfInitRealImag[0])), 
    1135             253 :                         nBandSize);
    1136             131 :             else if( !CPLIsNan(adfInitRealImag[0]) && adfInitRealImag[0] == 0.0 &&
    1137              43 :                      !CPLIsNan(adfInitRealImag[1]) && adfInitRealImag[1] == 0.0 )
    1138                 :             {
    1139              43 :                 memset( pBandData, 0, nBandSize );
    1140                 :             }
    1141               4 :             else if( !CPLIsNan(adfInitRealImag[1]) && adfInitRealImag[1] == 0.0 )
    1142                 :             {
    1143                 :                 GDALCopyWords( &adfInitRealImag, GDT_Float64, 0, 
    1144                 :                                pBandData,psWO->eWorkingDataType,nWordSize,
    1145               2 :                                nBlockXSize * nBlockYSize );
    1146                 :             }
    1147                 :             else
    1148                 :             {
    1149                 :                 GDALCopyWords( &adfInitRealImag, GDT_CFloat64, 0, 
    1150                 :                                pBandData,psWO->eWorkingDataType,nWordSize,
    1151               0 :                                nBlockXSize * nBlockYSize );
    1152                 :             }
    1153                 :         }
    1154                 : 
    1155              99 :         CSLDestroy( papszInitValues );
    1156                 :     }
    1157                 : 
    1158                 : /* -------------------------------------------------------------------- */
    1159                 : /*      Warp into this buffer.                                          */
    1160                 : /* -------------------------------------------------------------------- */
    1161                 :     CPLErr eErr;
    1162                 : 
    1163                 :     eErr = 
    1164                 :         poWarper->WarpRegionToBuffer( 
    1165                 :             iBlockX * nBlockXSize, iBlockY * nBlockYSize, 
    1166                 :             nBlockXSize, nBlockYSize,
    1167              99 :             pabyDstBuffer, psWO->eWorkingDataType );
    1168                 : 
    1169              99 :     if( eErr != CE_None )
    1170                 :     {
    1171               0 :         VSIFree( pabyDstBuffer );
    1172               0 :         return eErr;
    1173                 :     }
    1174                 :                         
    1175                 : /* -------------------------------------------------------------------- */
    1176                 : /*      Copy out into cache blocks for each band.                       */
    1177                 : /* -------------------------------------------------------------------- */
    1178             227 :     for( iBand = 0; iBand < MIN(nBands, psWO->nBandCount); iBand++ )
    1179                 :     {
    1180                 :         GDALRasterBand *poBand;
    1181                 :         GDALRasterBlock *poBlock;
    1182                 : 
    1183             128 :         poBand = GetRasterBand(iBand+1);
    1184             128 :         poBlock = poBand->GetLockedBlockRef( iBlockX, iBlockY, TRUE );
    1185                 : 
    1186             128 :         if( poBlock != NULL )
    1187                 :         {
    1188             128 :             if ( poBlock->GetDataRef() != NULL )
    1189                 :             {
    1190                 :                 GDALCopyWords( pabyDstBuffer + iBand*nBlockXSize*nBlockYSize*nWordSize,
    1191                 :                             psWO->eWorkingDataType, nWordSize, 
    1192                 :                             poBlock->GetDataRef(), 
    1193                 :                             poBlock->GetDataType(), 
    1194                 :                             GDALGetDataTypeSize(poBlock->GetDataType())/8,
    1195             128 :                             nBlockXSize * nBlockYSize );
    1196                 :             }
    1197                 : 
    1198             128 :             poBlock->DropLock();
    1199                 :         }
    1200                 :     }
    1201                 : 
    1202              99 :     VSIFree( pabyDstBuffer );
    1203                 :     
    1204              99 :     return CE_None;
    1205                 : }
    1206                 : 
    1207                 : /************************************************************************/
    1208                 : /*                              AddBand()                               */
    1209                 : /************************************************************************/
    1210                 : 
    1211              27 : CPLErr VRTWarpedDataset::AddBand( GDALDataType eType, char **papszOptions )
    1212                 : 
    1213                 : {
    1214                 :     UNREFERENCED_PARAM( papszOptions );
    1215                 : 
    1216                 :     SetBand( GetRasterCount() + 1,
    1217              27 :              new VRTWarpedRasterBand( this, GetRasterCount() + 1, eType ) );
    1218                 : 
    1219              27 :     return CE_None;
    1220                 : }
    1221                 : 
    1222                 : /************************************************************************/
    1223                 : /* ==================================================================== */
    1224                 : /*                        VRTWarpedRasterBand                           */
    1225                 : /* ==================================================================== */
    1226                 : /************************************************************************/
    1227                 : 
    1228                 : /************************************************************************/
    1229                 : /*                        VRTWarpedRasterBand()                         */
    1230                 : /************************************************************************/
    1231                 : 
    1232             110 : VRTWarpedRasterBand::VRTWarpedRasterBand( GDALDataset *poDS, int nBand,
    1233             110 :                                           GDALDataType eType )
    1234                 : 
    1235                 : {
    1236             110 :     Initialize( poDS->GetRasterXSize(), poDS->GetRasterYSize() );
    1237                 : 
    1238             110 :     this->poDS = poDS;
    1239             110 :     this->nBand = nBand;
    1240             110 :     this->eAccess = GA_Update;
    1241                 : 
    1242                 :     ((VRTWarpedDataset *) poDS)->GetBlockSize( &nBlockXSize, 
    1243             110 :                                                &nBlockYSize );
    1244                 : 
    1245             110 :     if( eType != GDT_Unknown )
    1246              29 :         this->eDataType = eType;
    1247             110 : }
    1248                 : 
    1249                 : /************************************************************************/
    1250                 : /*                        ~VRTWarpedRasterBand()                        */
    1251                 : /************************************************************************/
    1252                 : 
    1253             110 : VRTWarpedRasterBand::~VRTWarpedRasterBand()
    1254                 : 
    1255                 : {
    1256             110 :     FlushCache();
    1257             110 : }
    1258                 : 
    1259                 : /************************************************************************/
    1260                 : /*                             IReadBlock()                             */
    1261                 : /************************************************************************/
    1262                 : 
    1263              99 : CPLErr VRTWarpedRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
    1264                 :                                      void * pImage )
    1265                 : 
    1266                 : {
    1267                 :     CPLErr eErr;
    1268              99 :     VRTWarpedDataset *poWDS = (VRTWarpedDataset *) poDS;
    1269                 :     GDALRasterBlock *poBlock;
    1270                 : 
    1271              99 :     poBlock = GetLockedBlockRef( nBlockXOff, nBlockYOff, TRUE );
    1272              99 :     if( poBlock == NULL )
    1273               0 :         return CE_Failure;
    1274                 : 
    1275              99 :     eErr = poWDS->ProcessBlock( nBlockXOff, nBlockYOff );
    1276                 : 
    1277              99 :     if( eErr == CE_None && pImage != poBlock->GetDataRef() )
    1278                 :     {
    1279                 :         int nDataBytes;
    1280                 :         nDataBytes = (GDALGetDataTypeSize(poBlock->GetDataType()) / 8)
    1281               0 :             * poBlock->GetXSize() * poBlock->GetYSize();
    1282               0 :         memcpy( pImage, poBlock->GetDataRef(), nDataBytes );
    1283                 :     }
    1284                 : 
    1285              99 :     poBlock->DropLock();
    1286                 : 
    1287              99 :     return eErr;
    1288                 : }
    1289                 : 
    1290                 : /************************************************************************/
    1291                 : /*                            IWriteBlock()                             */
    1292                 : /************************************************************************/
    1293                 : 
    1294              16 : CPLErr VRTWarpedRasterBand::IWriteBlock( int nBlockXOff, int nBlockYOff,
    1295                 :                                      void * pImage )
    1296                 : 
    1297                 : {
    1298                 :     CPLErr eErr;
    1299              16 :     VRTWarpedDataset *poWDS = (VRTWarpedDataset *) poDS;
    1300                 : 
    1301                 :     /* This is a bit tricky. In the case we are warping a VRTWarpedDataset */
    1302                 :     /* with a destination alpha band, IWriteBlock can be called on that alpha */
    1303                 :     /* band by GDALWarpDstAlphaMasker */
    1304                 :     /* We don't need to do anything since the data will be kept in the block */
    1305                 :     /* cache by VRTWarpedRasterBand::IReadBlock */
    1306              16 :     if (poWDS->poWarper->GetOptions()->nDstAlphaBand == nBand)
    1307                 :     {
    1308              16 :         eErr = CE_None;
    1309                 :     }
    1310                 :     else
    1311                 :     {
    1312                 :         /* Otherwise, call the superclass method, that will fail of course */
    1313               0 :         eErr = VRTRasterBand::IWriteBlock(nBlockXOff, nBlockYOff, pImage);
    1314                 :     }
    1315                 : 
    1316              16 :     return eErr;
    1317                 : }
    1318                 : 
    1319                 : /************************************************************************/
    1320                 : /*                              XMLInit()                               */
    1321                 : /************************************************************************/
    1322                 : 
    1323              81 : CPLErr VRTWarpedRasterBand::XMLInit( CPLXMLNode * psTree, 
    1324                 :                                   const char *pszVRTPath )
    1325                 : 
    1326                 : {
    1327              81 :     return VRTRasterBand::XMLInit( psTree, pszVRTPath );
    1328                 : }
    1329                 : 
    1330                 : /************************************************************************/
    1331                 : /*                           SerializeToXML()                           */
    1332                 : /************************************************************************/
    1333                 : 
    1334              15 : CPLXMLNode *VRTWarpedRasterBand::SerializeToXML( const char *pszVRTPath )
    1335                 : 
    1336                 : {
    1337              15 :     CPLXMLNode *psTree = VRTRasterBand::SerializeToXML( pszVRTPath );
    1338                 : 
    1339                 : /* -------------------------------------------------------------------- */
    1340                 : /*      Set subclass.                                                   */
    1341                 : /* -------------------------------------------------------------------- */
    1342                 :     CPLCreateXMLNode( 
    1343                 :         CPLCreateXMLNode( psTree, CXT_Attribute, "subClass" ), 
    1344              15 :         CXT_Text, "VRTWarpedRasterBand" );
    1345                 : 
    1346              15 :     return psTree;
    1347                 : }
    1348                 : 
    1349                 : /************************************************************************/
    1350                 : /*                          GetOverviewCount()                          */
    1351                 : /************************************************************************/
    1352                 : 
    1353              11 : int VRTWarpedRasterBand::GetOverviewCount()
    1354                 : 
    1355                 : {
    1356              11 :     VRTWarpedDataset *poWDS = (VRTWarpedDataset *) poDS;
    1357                 : 
    1358              11 :     return poWDS->nOverviewCount;
    1359                 : }
    1360                 : 
    1361                 : /************************************************************************/
    1362                 : /*                            GetOverview()                             */
    1363                 : /************************************************************************/
    1364                 : 
    1365               1 : GDALRasterBand *VRTWarpedRasterBand::GetOverview( int iOverview )
    1366                 : 
    1367                 : {
    1368               1 :     VRTWarpedDataset *poWDS = (VRTWarpedDataset *) poDS;
    1369                 : 
    1370               1 :     if( iOverview < 0 || iOverview >= poWDS->nOverviewCount )
    1371               0 :         return NULL;
    1372                 :     else
    1373               1 :         return poWDS->papoOverviews[iOverview]->GetRasterBand( nBand );
    1374                 : }
    1375                 : 

Generated by: LCOV version 1.7