State-wide DRCM Processing

or: why big mosaics are cool

Problem:

• Wanted to know how many valid seasons per pixel in NT/Qld
• Too lazy to recall images, count pixels, mosaic and summarise
• Why not use the Seasonal Fractional Cover Mosaics?
• With rios, how hard can that be?

Plotted time vs file size.

Plotted time vs file size.

Plotted time vs file size.

Be Careful with Both creation options and rios options

    controls = applier.ApplierControls()
    controls.setOutputDriverName('GTiff')
    options = ['COMPRESS=LZW', 'BIGTIFF=YES', 'TILED=YES', 
               'INTERLEAVE=BAND', 
               'BLOCKXSIZE=256', 'BLOCKYSIZE=256']
    controls.setCreationOptions(options)                                      
    controls.setWindowXsize(256)
    controls.setWindowYsize(256)
    

median no. of seasons is 92 (66,105). Maximum possible is 112.

Writing in c

• faster (up to 10x)
• more difficult (but with mosaics not that difficult)
• All mosaics same size (usually), same grid
• seaonal ones masked and ready to go.

Useful for DRCM

• uses a huge (1501x1501) smoothing window
• edge effects potentially large
• complicates division and mosaicing of tiles.
• convenient if we can do in one go

    // assume all bands same size
    nXSize = GDALGetRasterBandXSize( fifthPercentileBand ); 
    nYSize = GDALGetRasterBandYSize( fifthPercentileBand ); 

    printf("allocation size: %lu\n",
          (unsigned long)nYSize* (unsigned long)nXSize*sizeof(uint8_t));
    
    // Allocate space
    pafScanline = (uint8_t *) CPLMalloc(nYSize*nXSize*sizeof(uint8_t));
    dixScanline = (uint8_t *) CPLMalloc(nYSize*nXSize*sizeof(uint8_t));
    refPixelScanline = (uint8_t *) CPLMalloc(nYSize*nXSize*sizeof(uint8_t));

    /* and read the whole thing in */
    GDALRasterIO(fifthPercentileBand , GF_Read, 0, 0, nXSize, nYSize,
		 pafScanline, nXSize, nYSize,  GDT_Byte,
		 0, 0 );
    GDALRasterIO(dixImageBand , GF_Read, 0, 0, nXSize, nYSize,
		 dixScanline, nXSize, nYSize,  GDT_Byte,
		 0, 0 );

    /* Create output */
    hDstDS = GDALCreate( hDriver, pszDstFilename, nXSize, nYSize, 1, GDT_Byte, 
                         papszOptions );
    destBand = GDALGetRasterBand( hDstDS, 1 );
    

Useful for DRCM

• requires estimate of 5th percentile of cover
• use 100 dix seasons (summer 1987 to sprint 2012)
• each dix about 3.7Gig in size
• process in c
• took 3hrs 45 minutes, used about 1.2gig of ram

fullsize

Useful for DRCM

• DRCM reference pixel takes longer
• have to move through image pixel by pixel
• ~ 15hours. Used about 8gig

fullsize

Useful for DRCM

• writing in c a possibiliy
  • quick routines still best in python
  • easier to maintain, read, fix
• but some situations c might be good:

We should forget about small efficiencies, say about 97% of the time: premature optimization is the root of all evil. Yet we should not pass up our opportunities in that critical 3%.

(Donald Knuth)