PostGIS  3.4.0dev-r@@SVN_REVISION@@

◆ RASTER_getPixelCentroids()

Datum RASTER_getPixelCentroids ( PG_FUNCTION_ARGS  )

Definition at line 664 of file rtpg_geometry.c.

665 {
666  FuncCallContext *funcctx;
667  TupleDesc tupdesc;
668  rt_pixel pix = NULL;
669  rt_pixel pix2;
670  int call_cntr;
671  int max_calls;
672  int i = 0;
673 
674  if (SRF_IS_FIRSTCALL()) {
675  MemoryContext oldcontext;
676  rt_pgraster *pgraster = NULL;
677  rt_raster raster = NULL;
678  rt_band band = NULL;
679  int nband = 1;
680  int numbands;
681  bool hasband = FALSE;
682  bool exclude_nodata_value = TRUE;
683  bool nocolumnx = TRUE;
684  bool norowy = TRUE;
685  int x = 0;
686  int y = 0;
687  int bounds[4] = {0};
688  int pixcount = 0;
689  double value = 0;
690  int isnodata = 0;
691 
692  LWPOINT* point = NULL;
693 
694  POSTGIS_RT_DEBUG(3, "RASTER_getPixelCentroids first call");
695 
696  /* create a function context for cross-call persistence */
697  funcctx = SRF_FIRSTCALL_INIT();
698 
699  /* switch to memory context appropriate for multiple function calls */
700  oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
701 
702  /* raster */
703  if (PG_ARGISNULL(0)) {
704  MemoryContextSwitchTo(oldcontext);
705  SRF_RETURN_DONE(funcctx);
706  }
707  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
708 
709  /* band */
710  if (!PG_ARGISNULL(1)) {
711  nband = PG_GETARG_INT32(1);
712  hasband = TRUE;
713  }
714 
715  /* column */
716  if (!PG_ARGISNULL(2)) {
717  bounds[0] = PG_GETARG_INT32(2);
718  bounds[1] = bounds[0];
719  nocolumnx = FALSE;
720  }
721 
722  /* row */
723  if (!PG_ARGISNULL(3)) {
724  bounds[2] = PG_GETARG_INT32(3);
725  bounds[3] = bounds[2];
726  norowy = FALSE;
727  }
728 
729  /* exclude NODATA */
730  if (!PG_ARGISNULL(4))
731  exclude_nodata_value = PG_GETARG_BOOL(4);
732 
733  /* deserialize raster */
734  raster = rt_raster_deserialize(pgraster, FALSE);
735  if (!raster) {
736  PG_FREE_IF_COPY(pgraster, 0);
737  ereport(ERROR, (
738  errcode(ERRCODE_OUT_OF_MEMORY),
739  errmsg("Could not deserialize raster")
740  ));
741  MemoryContextSwitchTo(oldcontext);
742  SRF_RETURN_DONE(funcctx);
743  }
744 
745  /* raster empty, return empty set */
746  if (rt_raster_is_empty(raster)) {
747  elog(NOTICE, "Raster is empty. Returning empty set");
748  rt_raster_destroy(raster);
749  PG_FREE_IF_COPY(pgraster, 0);
750  MemoryContextSwitchTo(oldcontext);
751  SRF_RETURN_DONE(funcctx);
752  }
753 
754  /* band specified, load band and info */
755  if (hasband) {
756  numbands = rt_raster_get_num_bands(raster);
757  POSTGIS_RT_DEBUGF(3, "band %d", nband);
758  POSTGIS_RT_DEBUGF(3, "# of bands %d", numbands);
759 
760  /* check band index */
761  if (nband < 1 || nband > numbands) {
762  elog(NOTICE, "Invalid band index (must use 1-based). Returning empty set");
763  rt_raster_destroy(raster);
764  PG_FREE_IF_COPY(pgraster, 0);
765  MemoryContextSwitchTo(oldcontext);
766  SRF_RETURN_DONE(funcctx);
767  }
768 
769  /* get band */
770  band = rt_raster_get_band(raster, nband - 1);
771  if (!band) {
772  elog(NOTICE, "Could not find band at index %d. Returning empty set", nband);
773 
774  rt_raster_destroy(raster);
775  PG_FREE_IF_COPY(pgraster, 0);
776 
777  MemoryContextSwitchTo(oldcontext);
778  SRF_RETURN_DONE(funcctx);
779  }
780 
781  if (!rt_band_get_hasnodata_flag(band))
782  exclude_nodata_value = FALSE;
783  }
784 
785  /* set bounds if columnx, rowy not set */
786  if (nocolumnx) {
787  bounds[0] = 1;
788  bounds[1] = rt_raster_get_width(raster);
789  }
790  if (norowy) {
791  bounds[2] = 1;
792  bounds[3] = rt_raster_get_height(raster);
793  }
794  POSTGIS_RT_DEBUGF(3, "bounds (min x, max x, min y, max y) = (%d, %d, %d, %d)",
795  bounds[0], bounds[1], bounds[2], bounds[3]);
796 
797  /* rowy */
798  pixcount = 0;
799  for (y = bounds[2]; y <= bounds[3]; y++) {
800  /* columnx */
801  for (x = bounds[0]; x <= bounds[1]; x++) {
802 
803  value = 0;
804  isnodata = TRUE;
805 
806  if (hasband) {
807  if (rt_band_get_pixel(band, x - 1, y - 1, &value, &isnodata) != ES_NONE) {
808  /* free already created centroid points */
809  for (i = 0; i < pixcount; i++)
810  lwgeom_free(pix[i].geom);
811  if (pixcount) pfree(pix);
812 
813  rt_band_destroy(band);
814  rt_raster_destroy(raster);
815  PG_FREE_IF_COPY(pgraster, 0);
816 
817  MemoryContextSwitchTo(oldcontext);
818  elog(ERROR, "RASTER_getPixelCentroids: Could not get pixel value");
819  SRF_RETURN_DONE(funcctx);
820  }
821 
822  /* don't continue if pixel is NODATA and to exclude NODATA */
823  if (isnodata && exclude_nodata_value) {
824  POSTGIS_RT_DEBUG(5, "pixel value is NODATA and exclude_nodata_value = TRUE");
825  continue;
826  }
827  }
828 
829  /* geometry */
830  point = rt_raster_pixel_as_centroid_point(raster, x - 1, y - 1);
831  if (!point) {
832  /* free already created centroid points */
833  for (i = 0; i < pixcount; i++)
834  lwgeom_free(pix[i].geom);
835  if (pixcount) pfree(pix);
836 
837  if (hasband) rt_band_destroy(band);
838  rt_raster_destroy(raster);
839  PG_FREE_IF_COPY(pgraster, 0);
840 
841  MemoryContextSwitchTo(oldcontext);
842  elog(ERROR, "RASTER_getPixelCentroids: Could not get pixel centroid");
843  SRF_RETURN_DONE(funcctx);
844  }
845 
846  /* allocate space for new point */
847  if (!pixcount) {
848  pix = palloc(sizeof(struct rt_pixel_t) * (pixcount + 1));
849  }
850  else {
851  pix = repalloc(pix, sizeof(struct rt_pixel_t) * (pixcount + 1));
852  }
853  if (pix == NULL) {
854  lwpoint_free(point);
855 
856  if (hasband) rt_band_destroy(band);
857  rt_raster_destroy(raster);
858  PG_FREE_IF_COPY(pgraster, 0);
859 
860  MemoryContextSwitchTo(oldcontext);
861  elog(ERROR, "RASTER_getPixelCentroids: Could not allocate memory for storing pixel centroids");
862  SRF_RETURN_DONE(funcctx);
863  }
864 
865  /* set pixel geometry */
866  pix[pixcount].geom = (LWGEOM *) point;
867  POSTGIS_RT_DEBUGF(5, "point @ %p", point);
868  POSTGIS_RT_DEBUGF(5, "geom @ %p", pix[pixcount].geom);
869 
870  /* set pixel coordinates */
871  pix[pixcount].x = x;
872  pix[pixcount].y = y;
873 
874  /* set pixel value */
875  pix[pixcount].value = value;
876 
877  /* NODATA */
878  if (hasband) {
879  if (exclude_nodata_value)
880  pix[pixcount].nodata = isnodata;
881  else
882  pix[pixcount].nodata = FALSE;
883  }
884  else {
885  pix[pixcount].nodata = isnodata;
886  }
887 
888  /* update pixcount*/
889  pixcount++;
890  }
891  }
892 
893  /* cleanup */
894  if (hasband) rt_band_destroy(band);
895  rt_raster_destroy(raster);
896  PG_FREE_IF_COPY(pgraster, 0);
897 
898  /* shortcut if no pixcount */
899  if (pixcount < 1) {
900  elog(NOTICE, "No pixels found for band %d", nband);
901  MemoryContextSwitchTo(oldcontext);
902  SRF_RETURN_DONE(funcctx);
903  }
904 
905  /* store pixel data */
906  funcctx->user_fctx = pix;
907 
908  /* set total number of tuples to be returned */
909  funcctx->max_calls = pixcount;
910  POSTGIS_RT_DEBUGF(3, "pixcount = %d", pixcount);
911 
912  /* build a tuple descriptor for our result type */
913  if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE) {
914  ereport(ERROR, (
915  errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
916  errmsg("function returning record called in context that cannot accept type record")
917  ));
918  }
919  /* set tuple descriptor */
920  BlessTupleDesc(tupdesc);
921  funcctx->tuple_desc = tupdesc;
922 
923  MemoryContextSwitchTo(oldcontext);
924  }
925 
926  /* stuff done on every call of the function */
927  funcctx = SRF_PERCALL_SETUP();
928 
929  call_cntr = funcctx->call_cntr;
930  max_calls = funcctx->max_calls;
931  tupdesc = funcctx->tuple_desc;
932  pix2 = funcctx->user_fctx;
933 
934  /* do when there is more left to send */
935  if (call_cntr < max_calls) {
936  Datum values[VALUES_LENGTH];
937  bool nulls[VALUES_LENGTH];
938  HeapTuple tuple;
939  Datum result;
940 
941  GSERIALIZED *gser = NULL;
942  size_t gser_size = 0;
943 
944  POSTGIS_RT_DEBUGF(3, "call number %d", call_cntr);
945 
946  memset(nulls, FALSE, sizeof(bool) * VALUES_LENGTH);
947 
948  /* convert LWGEOM go GSERIALIZED */
949  gser = gserialized_from_lwgeom(pix2[call_cntr].geom, &gser_size);
950  lwgeom_free(pix2[call_cntr].geom); /* free geometry object */
951 
952  /* set result tuple data */
953  values[0] = PointerGetDatum(gser);
954  if (pix2[call_cntr].nodata)
955  nulls[1] = TRUE;
956  else
957  values[1] = Float8GetDatum(pix2[call_cntr].value);
958  values[2] = Int32GetDatum(pix2[call_cntr].x);
959  values[3] = Int32GetDatum(pix2[call_cntr].y);
960 
961  /* build a tuple */
962  tuple = heap_form_tuple(tupdesc, values, nulls);
963 
964  /* make the tuple into a datum */
965  result = HeapTupleGetDatum(tuple);
966 
967  SRF_RETURN_NEXT(funcctx, result);
968  }
969  /* do when there is no more left */
970  else {
971  pfree(pix2);
972  SRF_RETURN_DONE(funcctx);
973  }
974 }
result
char result[OUT_DOUBLE_BUFFER_SIZE]
Definition: cu_print.c:262
TRUE
#define TRUE
Definition: dbfopen.c:73
FALSE
#define FALSE
Definition: dbfopen.c:72
gserialized_from_lwgeom
GSERIALIZED * gserialized_from_lwgeom(LWGEOM *geom, size_t *size)
Allocate a new GSERIALIZED from an LWGEOM.
Definition: gserialized.c:222
lwpoint_free
void lwpoint_free(LWPOINT *pt)
Definition: lwpoint.c:213
lwgeom_free
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1155
rt_band_get_hasnodata_flag
int rt_band_get_hasnodata_flag(rt_band band)
Get hasnodata flag value.
Definition: rt_band.c:674
rt_band_get_pixel
rt_errorstate rt_band_get_pixel(rt_band band, int x, int y, double *value, int *nodata)
Get pixel value.
Definition: rt_band.c:1376
rt_raster_destroy
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:86
rt_raster_pixel_as_centroid_point
LWPOINT * rt_raster_pixel_as_centroid_point(rt_raster rast, int x, int y)
Get a raster pixel centroid point.
Definition: rt_geometry.c:668
ES_NONE
@ ES_NONE
Definition: librtcore.h:182
rt_band_destroy
void rt_band_destroy(rt_band band)
Destroy a raster band.
Definition: rt_band.c:340
rt_raster_get_num_bands
uint16_t rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:376
rt_raster_get_height
uint16_t rt_raster_get_height(rt_raster raster)
Definition: rt_raster.c:133
rt_raster_get_width
uint16_t rt_raster_get_width(rt_raster raster)
Definition: rt_raster.c:125
rt_raster_deserialize
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_serialize.c:725
rt_raster_is_empty
int rt_raster_is_empty(rt_raster raster)
Return TRUE if the raster is empty.
Definition: rt_raster.c:1362
rt_raster_get_band
rt_band rt_raster_get_band(rt_raster raster, int bandNum)
Return Nth band, or NULL if unavailable.
Definition: rt_raster.c:385
genraster.value
int value
Definition: genraster.py:62
ovdump.band
band
Definition: ovdump.py:58
pixval.nband
nband
Definition: pixval.py:53
rtrowdump.raster
raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:121
VALUES_LENGTH
#define VALUES_LENGTH
Definition: rtpg_geometry.c:658
POSTGIS_RT_DEBUG
#define POSTGIS_RT_DEBUG(level, msg)
Definition: rtpostgis.h:65
POSTGIS_RT_DEBUGF
#define POSTGIS_RT_DEBUGF(level, msg,...)
Definition: rtpostgis.h:69
rt_band_t
Definition: librtcore.h:2466
rt_pixel_t::y
int y
Definition: librtcore.h:2488
rt_pixel_t::geom
LWGEOM * geom
Definition: librtcore.h:2493
rt_pixel_t::value
double value
Definition: librtcore.h:2491
rt_pixel_t::nodata
uint8_t nodata
Definition: librtcore.h:2490
rt_pixel_t::x
int x
Definition: librtcore.h:2487
rt_pixel_t
Definition: librtcore.h:2486
rt_raster_serialized_t
Struct definitions.
Definition: librtcore.h:2403
rt_raster_t
Definition: librtcore.h:2437

References ovdump::band, ES_NONE, FALSE, rt_pixel_t::geom, gserialized_from_lwgeom(), lwgeom_free(), lwpoint_free(), pixval::nband, rt_pixel_t::nodata, POSTGIS_RT_DEBUG, POSTGIS_RT_DEBUGF, rtrowdump::raster, result, rt_band_destroy(), rt_band_get_hasnodata_flag(), rt_band_get_pixel(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_get_height(), rt_raster_get_num_bands(), rt_raster_get_width(), rt_raster_is_empty(), rt_raster_pixel_as_centroid_point(), TRUE, rt_pixel_t::value, genraster::value, VALUES_LENGTH, rt_pixel_t::x, pixval::x, rt_pixel_t::y, and pixval::y.

Here is the call graph for this function: