/*
 * Separate calculating projected coordinates and transforming the pixels.
 *
 * Copyright (C) 2018  Antonio Ospite <ao2@ao2.it>
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * The idea behind these functions is that in some use cases the perspective
 * deformation does not change, it's only the projected data which changes, so
 * there's no need to recalculate the projection parameters every time.
 *
 * It is enough to apply the pre-calculated projection to the new data,
 * possibly using a region-of-interest, to avoid clipping the image in case
 * it's larger than the viewport.
 */

#include "projective_split.h"

/*
 * Pre-calculate the transformation map.
 *
 * This can be useful if the transform mapping does not change but we want to
 * apply the same transform to different data.
 */
struct p *
_pixProjectiveSampled_precalc_map(BOX        *viewport,
					l_float32  *vc)
{
	l_int32    i, j, x, y;
	struct p   *map;

	if (!vc)
		return (struct p *)ERROR_PTR("vc not defined", __func__, NULL);

	map = malloc(viewport->w * viewport->h * sizeof(*map));
	if (map == NULL)
		return (struct p *)ERROR_PTR("cannot allocate map", __func__, NULL);

	for (i = 0; i < viewport->h; i++) {
		for (j = 0; j < viewport->w; j++) {
			/* XXX: transform i and j by considering the viewport origin */
			projectiveXformSampledPt(vc, j, i, &x, &y);

			/* XXX: This is a hack to avoid conditionals when
			 * applying the map.
			 * Basically we set all out of bound pixels to be the same
			 * as the pixel at 0,0.
			 *
			 * When applying the map make sure that the pixel at
			 * 0,0 is set to the desired "background" color.
			 */
			if (x < 0 || y < 0 || x >= viewport->w || y >= viewport->h) {
				x = 0;
				y = 0;
			}

			map[i * viewport->w + j].x = x;
			map[i * viewport->w + j].y = y;
		}
	}

	return map;
}

#ifdef USE_NEON

#include <arm_neon.h>

/*
 * Apply the map pre-calculated transformation map to an image, but only to
 * a region of interest.
 *
 * This is a NEON optimized version, however the speedup is not that
 * impressive because moving pixels around always means accessing uncontiguous
 * memory.
 */
PIX *
_pixProjectiveSampled_apply_map_dest_roi(PIX        *pixs,
			       struct p   *map,
			       l_int32     incolor,
			       PIX        *pixd,
			       BOX        *roi)
{
	int32_t     i, j, w, h, d, wpls, wpld;
	uint32_t   *datas, *datad, *lined;
	l_uint32 pixel0;

	uint32_t stride;
	uint32x4x2_t point;
	uint32x4_t lines;
	uint32x4_t roi_x;
	uint32x4_t roi_y;
	uint32x4_t x_offset;
	uint32x4_t y_offset;
	uint32x4_t wpls_v;

	PROCNAME("pixProjectiveSampled_neon");

	if (!pixs)
		return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
	if (!map)
		return (PIX *)ERROR_PTR("map not defined", procName, NULL);
	if (incolor != L_BRING_IN_WHITE && incolor != L_BRING_IN_BLACK)
		return (PIX *)ERROR_PTR("invalid incolor", procName, NULL);
	pixGetDimensions(pixs, &w, &h, &d);
	if (d != 32)
		return (PIX *)ERROR_PTR("depth not 32", procName, NULL);
	if (!pixd)
		return (PIX *)ERROR_PTR("pixd not defined", procName, NULL);

	/* Scan over the dest pixels */
	datas = pixGetData(pixs);
	wpls = pixGetWpl(pixs);
	datad = pixGetData(pixd);
	wpld = pixGetWpl(pixd);

	roi_x = vmovq_n_u32(roi->x);
	roi_y = vmovq_n_u32(roi->y);
	wpls_v = vmovq_n_u32(wpls);

	/*
	 * Save the value of the pixel at 0,0 in the destination image,
	 * and set it to black, because the map uses the pixel at 0,0 for out
	 * of bounds pixels.
	 */
	pixGetPixel(pixs, roi->x, roi->y, &pixel0);
	pixSetPixel(pixs, roi->x, roi->y, 0);

	for (i = 0; i < roi->h; i++) {
		lined = datad + i * wpld;
		stride = i * roi->w;
		for (j = 0; j < roi->w; j += 4) {
			unsigned int map_index = stride + j;

			point = vld2q_u32((uint32_t *)(map + map_index));

			x_offset = vaddq_u32(point.val[0], roi_x);

			y_offset = vaddq_u32(point.val[1], roi_y);

			lines = vmlaq_u32(x_offset, y_offset, wpls_v);

			lined[j + 0] = *(datas + lines[0]);
			lined[j + 1] = *(datas + lines[1]);
			lined[j + 2] = *(datas + lines[2]);
			lined[j + 3] = *(datas + lines[3]);
		}
	}

	/* restore the previous value */
	pixSetPixel(pixs, roi->x, roi->y, pixel0);

	return pixd;
}

#else

/*
 * Apply the map pre-calculated transformation map to an image, but only to
 * a region of interest.
 */
PIX *
_pixProjectiveSampled_apply_map_dest_roi(PIX        *pixs,
			       struct p   *map,
			       l_int32     incolor,
			       PIX        *pixd,
			       BOX        *roi)
{
	l_int32     i, j, w, h, d, wpls, wpld;
	l_uint32   *datas, *datad, *lined;
	l_uint32 datas_offset;
	l_uint32 stride;
	l_uint32 pixel0;
	struct p *point;

	PROCNAME("pixProjectiveSampled");

	if (!pixs)
		return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
	if (!map)
		return (PIX *)ERROR_PTR("map not defined", procName, NULL);
	if (incolor != L_BRING_IN_WHITE && incolor != L_BRING_IN_BLACK)
		return (PIX *)ERROR_PTR("invalid incolor", procName, NULL);

	pixGetDimensions(pixs, &w, &h, &d);
	if (d != 32)
		return (PIX *)ERROR_PTR("depth not 32", procName, NULL);
	if (!pixd)
		return (PIX *)ERROR_PTR("pixd not defined", procName, NULL);

	/* Scan over the dest pixels */
	datas = pixGetData(pixs);
	wpls = pixGetWpl(pixs);
	datad = pixGetData(pixd);
	wpld = pixGetWpl(pixd);

	/*
	 * Save the value of the pixel at 0,0 in the destination image,
	 * and set it to black, because the map uses the pixel at 0,0 for out
	 * of bounds pixels.
	 *
	 * This is just a dirty trick to avoid a conditional deep in the loop
	 * below.
	 *
	 * A more general algorithm would check for out of bounds explicitly,
	 * with something like:
	 *
	 *	if (point->x >= 0 && point->y >= 0 && point->x < roi->w && point->y < roi->h) {
	 *		datas_offset = (point->x + roi->x) + (point->y + roi->y) * wpls;
	 *		*(lined + j) = *(datas + datas_offset);
	 *	}
	 */
	pixGetPixel(pixs, roi->x, roi->y, &pixel0);
	pixSetPixel(pixs, roi->x, roi->y, 0);

	for (i = 0; i < roi->h; i++) {
		lined = datad + i * wpld;
		stride = i * roi->w;
		for (j = 0; j < roi->w; j++) {
			point = map + stride + j;

			datas_offset = (point->x + roi->x) + (point->y + roi->y) * wpls;
			*(lined + j) = *(datas + datas_offset);
		}
	}

	/* restore the previous value */
	pixSetPixel(pixs, roi->x, roi->y, pixel0);

	return pixd;
}


#endif