2

I have a list of rasters some with the same date that I have to turn into VRTs with the same date too. The problem is that when I create a loop in the Gdal list it overwrites the VRTs that have the same name. While I want it to add them in the same VRT. E.g:

root = 'F:...\\MSAVI'
out_vrt = 'F:...\\out_vrt\\'

raster = [os.path.join(root, file) for root, dirs, files in os.walk(root) 
               for file in files if os.path.splitext(file)[1] == '.tif' ]
raster 
['.../33TVE6/S2_20190712_079_33TVE6_B.tif',
 '.../32TQQ3/S2_20190712_079_33TQQ3_B.tif',
  ...] 

dates = []
for file in raster:
    dates.append(file.split('\\')[6].split('_')[1])

dates
['20190817',
 '20190901',
 ...]


for day in dates:
    for file in raster:
        data_path = file.split('\\')[6].split('_')[1]
        name = 'S2x2A_' + day + '_MSAVI.vrt'
        tile = file.split('\\')[6].split('_')[3]
        
        if day in data_path:

            my_vrt = gdal.BuildVRT(out_vrt+ name, file)
            my_vrt = None

The output is a VRT with only one raster:

<VRTDataset rasterXSize="4001" rasterYSize="4001">
  <SRS dataAxisToSRSAxisMapping="1,2">PROJCS["WGS 84 / UTM zone 33N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",15],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32633"]]</SRS>
  <GeoTransform>  4.3500000000000000e+05,  1.0000000000000000e+01,  0.0000000000000000e+00,  4.4983300000000000e+06,  0.0000000000000000e+00, -1.0000000000000000e+01</GeoTransform>
  <VRTRasterBand dataType="Int16" band="1">
    <NoDataValue>-32768</NoDataValue>
    <ColorInterp>Gray</ColorInterp>
    <ComplexSource>
      <SourceFilename relativeToVRT="0">F:...\33TVE6\S2_20190712_079_33TVE6_B_MSAVI.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="4001" RasterYSize="4001" DataType="Int16" BlockXSize="4001" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="4001" ySize="4001" />
      <DstRect xOff="0" yOff="0" xSize="4001" ySize="4001" />
      <NODATA>-32768</NODATA>
    </ComplexSource>
  </VRTRasterBand>
</VRTDataset>

How can I add rasters with the same date to get a VRT with more rasters without overwriting the previous VRT?

3

BuildVRT can take a list of files as it's 2nd (source) parameter so you need to find all the files that match your date (using a full path) and pass that list to the function.

So something like:

for day in dates:
  f = [file for file in rasters if rasters.contains(day)]
  my_vrt = gdal.BuildVRT(out_vrt+name, f)
| improve this answer | |
  • It works, good. If I wanted to divide into 'f' the rasters that differ for the '33TVE6' or '32TQQ3' tile and get a list with only files from spindle 32 or spindle 33, how could I do it? – vins_26 Aug 23 at 15:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.