def write_blockwise(input_raster, output_raster):
    # input_raster and output_raster are gdal class object
    """
    input_raster = gdal.Open(path_input)
    output_raster = gdal.Open(path_output)
    """
    # Geting dimension of raster
    dimension = input_raster.num_band

    # Iterating over number of bands
    for num_band in range(dimension):
        print('Processing: Band %d' % (num_band))
        input_band = input_raster.ds.GetRasterBand(num_band+1)
        output_band = output_raster.GetRasterBand(num_band+1)
        block_x, block_y = input_band.GetBlockSize()

        # We can increase blocksize according our suitability also. It makes no difference
        # block_x, block_y = 5000, 5000

         # Step size
        size_x = input_band.XSize
        size_y = input_band.YSize

        # Iterating over total blocks with blocksize as step
        for x in range(0, int(size_x), int(block_x)):
            if x + block_x < size_x:
                col = block_x
            else:
                col = size_x - x

            for y in range(0, int(size_y), int(block_y)):
                if y + block_y < size_y:
                    row = block_y
                else:
                    row = size_y - y

                # Reading specific raster from large TIF
                array = input_band.ReadAsArray(x, y, col, row)

                """
                  To Do some fancy analysis
                """
                output_band.WriteArray(array, x, y)
                array = None