[gdal-dev] Amplitude virtual bands for complex datasets ?

Wed Jul 6 00:15:13 PDT 2016

Hi all,

I hacked the SetBand() method so as to report the new medata in 
DERIVED_SUBDATASETS mdd. Now if I gdal info a S1 dataset I get :

$ gdalinfo -nogcp -mdd DERIVED_SUBDATASETS 
s1a-s6-slc-vv-20150619t195043-20150619t195101-006447-00887d-001.tiff
Driver: GTiff/GeoTIFF
Files: s1a-s6-slc-vv-20150619t195043-20150619t195101-006447-00887d-001.tiff
s1a-s6-slc-vv-20150619t195043-20150619t195101-006447-00887d-001.tiff.aux.xml
Size is 17663, 31106
[...]
Metadata (DERIVED_SUBDATASETS):
   DERIVED_SUBDATASET_1_DESC=Complex amplitude of bands from 
/media/michelj/HUGE/LargeInput/SENTINEL1/S1A_S6_SLC__1SSV_20150619T195043/measurement/s1a-s6-slc-vv-20150619t195043-20150619t195101-006447-00887d-001.tiff
DERIVED_SUBDATASET_1_NAME=DERIVED_SUBDATASET:COMPLEX_AMPLITUDE:/media/michelj/HUGE/LargeInput/SENTINEL1/S1A_S6_SLC__1SSV_20150619T195043/measurement/s1a-s6-slc-vv-20150619t195043-20150619t195101-006447-00887d-001.tiff
[...]
Band 1 Block=17663x1 Type=CInt16, ColorInterp=Gray

Even, I think the AddBand() method is obviously not the right place to 
do that, but I do not get how tweak the GetMetadata() / 
GetMetadataDomainList(). As far as I can tell from the code, datasets 
never directly declare new domains by hacking GetMetadataDomainList(). 
Calling SetMetadataItem() with a new domain is enough to create it (that 
is what I did). But I surely miss something here.

Also, can we discuss syntax (COMPLEX_AMPLITUDE is not very meaningful) 
and behaviour (for now all bands are considered even if they are not 
complex).

Last, it does not work for real subdatasets. For instance, for S2 data :

$ gdalinfo -mdd DERIVED_SUBDATASETS 
SENTINEL2_L1C:S2A_OPER_MTD_SAFL1C_PDMC_20150519T113602_R070_V20130707T172156_20130707T172156.xml:20m:EPSG_32615

Driver: SENTINEL2/Sentinel 2
[...]
Metadata (DERIVED_SUBDATASETS):
   DERIVED_SUBDATASET_1_DESC=Complex amplitude of bands from
   DERIVED_SUBDATASET_1_NAME=DERIVED_SUBDATASET:COMPLEX_AMPLITUDE:
[...]

But calling
$ gdalinfo 
DERIVED_SUBDATASET:COMPLEX_AMPLITUDE:SENTINEL2_L1C:S2A_OPER_MTD_SAFL1C_PDMC_20150519T113602_R070_V20130707T172156_20130707T172156.xml:20m:EPSG_32615

seems to work.

I will try to do the pixel function proper integration on a separate 
branch, so as to avoid mixing different purposes.

Regards,

Julien

Le 04/07/2016 à 15:04, Even Rouault a écrit :
> Le lundi 04 juillet 2016 14:30:26, Julien Michel a écrit :
>> Hi,
>>
>> Started the work here :
>>
>> https://github.com/OSGeo/gdal/compare/trunk...jmichel-otb:enhance-complex-d
>> atasets
>>
>> For now, recognizes the following syntax :
>>
>> gdal_translate -srcwin 1000 1000 1000 1000
>> DERIVED_SUBDATASET:COMPLEX_AMPLITUDE:s1a-s6-slc-vv-20150619t195043-20150619
>> t195101-006447-00887d-001.tiff ~/tmp/test.tif
>>
>> I need to see about exposing those derived datasets properly.
>>
>> Antonio, I copied one of your pixel function for fast prototyping, I
>> hope you do not mind. Of course upon completion of this work we should
>> merge properly all your pixel functions with proper
>> credits/copyright/licence. Thing is I did not know where to put them
>> within gdal.
> frmts/vrt seem to be the appropriate place, and make the VRTDriver register
> them.
>
>> Regards,
>>
>> Julien
>>
>> Le 04/07/2016 à 11:06, Even Rouault a écrit :
>>> Hi,
>>>
>>>> To check if I understood well :
>>>>
>>>> I will create a driver that will recognize the
>>>> "DERIVED_SUBDATASET:Amplitude:original_datasetname" syntax. This driver
>>>> needs to now if "original_datasetname" is of complex type (to report it
>>>> can open it)
>>> That or we could also possibly accept non-complex types (assume imaginary
>>> part = 0). Antonio's functions handle complex & non-complex AFAICS
>>>
>>>> , and will also expose the derived band somehow (can I
>>>> benefit from the CreateDerivedBand API,
>>> If you write it, yes... (it doesn't exist yet if I wasn't clear)
>>>
>>> If you follow the VRT pixel function road, that's a matter of building a
>>> on- the-fly VRT dataset.
>>>
>>> http://www.gdal.org/gdal_vrttut.html#gdal_vrttut_creation could possibly
>>> be used. But AFAICS setting the pixel function name isn't supported.
>>> Could be done by implementing CPLErr
>>> VRTSourcedRasterBand::SetMetadataItem( const char *pszName, const char
>>> *pszValue, const char *pszDomain ), similarly to what
>>> VRTSourcedRasterBand::SetMetadataItem(...) handles for vrt sources.
>>>
>>>> and the pixel functions from Antonio ?).
>>> Would make sense to internalize this code (uses X/MIT license)
>>>
>>>> Also, "original_datasetname" can be a subdataset itself.
>>> Yes, that shouldn't matter to the driver. It is just a string that it
>>> provides to GDALOpen() (some care must be taken when parsing
>>> DERIVED_SUBDATASET:algname:original_datasetname, since
>>> original_datasetname can have columns too)
>>>
>>>> There should also be a mechanism that will report the existing
>>>> DERIVED_SUBDATASETS upon query (this I did not get yet).
>>> Implement in GDALDataset, the following method (from GDALMajorObject)
>>>
>>> virtual char **GetMetadata( const char * pszDomain = "" );
>>>
>>> char **GDALDataset::GetMetadata( const char * pszDomain = "" )
>>> {
>>>
>>> 	if( pszDomain != NULL && EQUAL(pszDomain, "DERIVED_SUBDATASET") )
>>> 	{
>>> 	
>>> 		// some stuff. Be careful that ownership of returned char** belongs
>>> 		// to the dataset object
>>> 	
>>> 	}
>>> 	else
>>> 	
>>> 		return GDALMajorObject::GetMetadata(pszDomain);
>>>
>>> }
>>>
>>> As well as :
>>>
>>> virtual char **GetMetadataDomainList(); to report DERIVED_SUBDATASET when
>>> it makes sense (contrary to GetMetadata(), ownership of the returned
>>> char** belongs to the caller)
>>>
>>> Even

-- 
Julien MICHEL
CNES - DCT/SI/AP - BPI 1219
18, avenue Edouard Belin
31401 Toulouse Cedex 09 - France
Tel: +33 561 282 894 - Fax: +33 561 283 109