tile_splitter.go

package main

import (
	"fmt"
	"image"
	"io/ioutil"
	"log"
	"os"
	"unsafe"

	"github.com/golang/snappy"
	"github.com/terrascope/geometry"
	"github.com/terrascope/proj4go"
	"github.com/terrascope/raster"
	"github.com/terrascope/scimage"
	"github.com/terrascope/scimage/cog"
	"github.com/terrascope/scimage/scicolor"
)

const (
	rasterName = "/Users/pablo/Downloads/tiles/cmrset_h%02dv%02d.tif"
	tileName   = "/Users/pablo/Downloads/tiles/cmrset_h%02dv%02ds%02d.snp"
	sinuProj   = "+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs "
	xSize      = 1111950.519666
	ySize      = 1111950.519667
	rasterSize = 2400
	tileSize   = 400
)

func getBlankRaster(h, v int) raster.Raster {
	x0 := float64(h-18) * xSize
	x1 := float64(h-17) * xSize
	y0 := float64(9-v) * ySize
	y1 := float64(8-v) * ySize

	return raster.Raster{scimage.NewBlankImage(scicolor.GrayS16Model{}, image.Rect(0, 0, 2400, 2400)),
		proj4go.Coverage{Proj4: sinuProj, BoundingBox: geometry.BBox(x0, y0, x1, y1)}}
}

func ReadRaster(h, v int) (raster.Raster, error) {
	fileName := fmt.Sprintf(rasterName, h, v)
	rast := getBlankRaster(h, v)

	rc, err := os.Open(fileName)
	if err != nil {
		return raster.Raster{}, fmt.Errorf("Error opening file: %s: %v", fileName, err)
	}

	rast.Image, err = cog.Decode(rc)
	rc.Close()
	if err != nil {
		return raster.Raster{}, fmt.Errorf("Error decoding image: %v", err)
	}

	return rast, nil
}

func GenerateTiles(h, v int, img image.Image) error {
	for i := 0; i < rasterSize/tileSize; i++ {
		for j := 0; j < rasterSize/tileSize; j++ {
			s := j*10 + i
			tName := fmt.Sprintf(tileName, h, v, s)

			rect := image.Rect(i*tileSize, j*tileSize, (i+1)*tileSize, (j+1)*tileSize)
			fmt.Println(rect)

			tile := img.(*scimage.GrayS16).SubImage(rect).(*scimage.GrayS16)
			pix := make([]int16, tileSize*tileSize)

			for ti := tile.Bounds().Min.X; ti < tile.Bounds().Max.X; ti++ {
				for tj := tile.Bounds().Min.Y; tj < tile.Bounds().Max.Y; tj++ {
					pix[(ti-tile.Bounds().Min.X)+(tj-tile.Bounds().Min.Y)*tileSize] = tile.GrayS16At(ti, tj).Y
				}
			}

			data := *(*[]byte)(unsafe.Pointer(&pix))
			fmt.Println(len(data))
			err := ioutil.WriteFile(tName, snappy.Encode(nil, data), 0644)
			if err != nil {
				return err
			}
		}
	}

	return nil
}

func main() {
	rast, err := ReadRaster(30, 12)
	if err != nil {
		log.Fatal(err)
	}

	im := rast.Image.(*scimage.GrayS16)
	GenerateTiles(30, 12, im)
}