Yet another intermediate release hudson-beamb-3-SUCCESS
authorAnders Henja <anders@baltrad.eu>
Wed, 23 Nov 2011 06:37:15 +0000 (07:37 +0100)
committerAnders Henja <anders@baltrad.eu>
Wed, 23 Nov 2011 06:37:15 +0000 (07:37 +0100)
lib/bbtopography.c
lib/bbtopography.h
lib/beamblockage.c
lib/beamblockagemap.c
lib/beamblockagemap.h

index 97c1d41..d32862a 100644 (file)
@@ -284,50 +284,6 @@ BBTopography_t* BBTopography_concatX(BBTopography_t* self, BBTopography_t* other
   return result;
 }
 
-int BBTopography_getLonLatFields(BBTopography_t* self, PolarScan_t* scan, RaveData2D_t** lon, RaveData2D_t** lat)
-{
-  int result = 0;
-  long nrays = 0, nbins = 0;
-  long ri = 0, bi = 0;
-  double elangle = 0.0;
-  RaveData2D_t *olon = NULL, *olat = NULL;
-
-  RAVE_ASSERT((self != NULL), "self == NULL");
-  RAVE_ASSERT((lon != NULL), "lon == NULL");
-  RAVE_ASSERT((lat != NULL), "lat == NULL");
-
-  if (scan == NULL) {
-    goto done;
-  }
-
-  nrays = PolarScan_getNrays(scan);
-  nbins = PolarScan_getNbins(scan);
-  elangle = PolarScan_getElangle(scan);
-
-  olon = RAVE_OBJECT_NEW(&RaveData2D_TYPE);
-  olat = RAVE_OBJECT_NEW(&RaveData2D_TYPE);
-
-  if (olon == NULL || olat == NULL) {
-    goto done;
-  }
-
-  if (!RaveData2D_createData(olon, nbins, nrays, RaveDataType_DOUBLE) ||
-      !RaveData2D_createData(olat, nbins, nrays, RaveDataType_DOUBLE)) {
-    goto done;
-  }
-
-  *lat = RAVE_OBJECT_COPY(olat);
-  *lon = RAVE_OBJECT_COPY(olon);
-
-
-  result = 1;
-done:
-  RAVE_OBJECT_RELEASE(olon);
-  RAVE_OBJECT_RELEASE(olat);
-
-  return result;
-}
-
 /*@} End of Interface functions */
 
 RaveCoreObjectType BBTopography_TYPE = {
index 2b1f472..c1ef68f 100644 (file)
@@ -204,6 +204,4 @@ int BBTopography_setValue(BBTopography_t* self, long col, long row, double value
  */
 BBTopography_t* BBTopography_concatX(BBTopography_t* self, BBTopography_t* other);
 
-int BBTopography_getLonLatFields(BBTopography_t* self, PolarScan_t* scan, RaveData2D_t** lon, RaveData2D_t** lat);
-
 #endif /* BBTOPOGRAPHY_H */
index 48e41df..58d1b6b 100644 (file)
@@ -148,20 +148,53 @@ const char* BeamBlockage_getTopo30Directory(BeamBlockage_t* self)
 RaveField_t* BeamBlockage_getBlockage(BeamBlockage_t* self, PolarScan_t* scan)
 {
   RaveField_t *field = NULL, *result = NULL;
-  double* ground_range = NULL;
+  RaveData2D_t *lon = NULL, *lat = NULL;
+  BBTopography_t* topo = NULL;
+  double lon_min = 0.0, lon_max = 0.0, lat_min = 0.0, lat_max = 0.0;
+
+  long ulx = 0, uly = 0, llx = 0, lly = 0;
 
   RAVE_ASSERT((self != NULL), "self == NULL");
   if (scan == NULL) {
     return NULL;
   }
-  ground_range = BeamBlockageInternal_computeGroundRange(self, scan);
-  if (ground_range == NULL) {
+
+  topo = BeamBlockageMap_readTopography(self->mapper,
+                                        PolarScan_getLatitude(scan),
+                                        PolarScan_getLongitude(scan),
+                                        PolarScan_getMaxDistance(scan));
+  if (topo == NULL) {
     goto done;
   }
 
+#ifdef KALLE
+  if (!BeamBlockageInternal_getLonLatFields(scan, &lon, &lat)) {
+    goto done;
+  }
+
+  if (!BeamBlockageInternal_min(lon, &lon_min) ||
+      !BeamBlockageInternal_max(lon, &lon_max) ||
+      !BeamBlockageInternal_min(lat, &lat_min) ||
+      !BeamBlockageInternal_max(lon, &lat_max)) {
+    goto done;
+  }
+  lon_min = floor((lon_min - BBTopography_getUlxmap(topo)) / BBTopography_getXDim(topo)) - 1;
+  lon_max = ceil((lon_max - BBTopography_getUlxmap(topo)) / BBTopography_getXDim(topo)) + 1;
+  lat_max = floor((BBTopography_getUlymap(topo) - lat_max) / BBTopography_getYDim(topo)) - 1;
+  lat_min = ceil((BBTopography_getUlymap(topo) - lat_min) / BBTopography_getYDim(topo)) + 1;
+
+#endif
+  /*
+  lon_min = floor( (min_double(lon,ri*ai)-ulxmap) / xdim ) - 1;
+  lon_max = ceil( (max_double(lon,ri*ai)-ulxmap) / xdim ) + 1;
+  lat_max = floor( (ulymap-max_double(lat,ri*ai)) / ydim ) - 1;
+  lat_min = ceil( (ulymap-min_double(lat,ri*ai)) / ydim ) + 1;
+  */
   result = RAVE_OBJECT_COPY(field);
 done:
   RAVE_OBJECT_RELEASE(field);
+  RAVE_OBJECT_RELEASE(lon);
+  RAVE_OBJECT_RELEASE(lat);
   return result;
 }
 
index ca32149..01e9d7f 100644 (file)
@@ -265,6 +265,172 @@ done:
   RAVE_OBJECT_RELEASE(field);
   return result;
 }
+/**
+ * Locates the min value in the rave data2d field.
+ * @param[in] field - the data field
+ * @param[out] omin - the minimum value
+ * @return 1 on success otherwise 0
+ */
+static int BeamBlockageMapInternal_min(RaveData2D_t* field, double* omin)
+{
+  double minv = 0.0, v = 0.0;
+  long xsize = 0, ysize = 0;
+  long x = 0, y = 0;
+  int result = 0;
+
+  RAVE_ASSERT((field != NULL), "field == NULL");
+  RAVE_ASSERT((omin != NULL), "omin == NULL");
+
+  xsize = RaveData2D_getXsize(field);
+  ysize = RaveData2D_getYsize(field);
+  if (x > 0 && y > 0) {
+    result = RaveData2D_getValue(field, 0, 0, &minv);
+    for (x = 0; result && x < xsize; x++) {
+      for (y = 0; result && y < ysize; y++) {
+        result = RaveData2D_getValue(field, x, y, &v);
+        if (v <= minv) {
+          minv = v;
+        }
+      }
+    }
+  }
+  *omin = minv;
+  return result;
+}
+
+/**
+ * Locates the max value in the rave data2d field.
+ * @param[in] field - the data field
+ * @param[out] omax - the minimum value
+ * @return 1 on success otherwise 0
+ */
+static int BeamBlockageMapInternal_max(RaveData2D_t* field, double* omax)
+{
+  double maxv = 0.0, v = 0.0;
+  long xsize = 0, ysize = 0;
+  long x = 0, y = 0;
+  int result = 0;
+
+  RAVE_ASSERT((field != NULL), "field == NULL");
+  RAVE_ASSERT((omax != NULL), "omax == NULL");
+
+  xsize = RaveData2D_getXsize(field);
+  ysize = RaveData2D_getYsize(field);
+  if (x > 0 && y > 0) {
+    result = RaveData2D_getValue(field, 0, 0, &maxv);
+    for (x = 0; result && x < xsize; x++) {
+      for (y = 0; result && y < ysize; y++) {
+        result = RaveData2D_getValue(field, x, y, &v);
+        if (v >= maxv) {
+          maxv = v;
+        }
+      }
+    }
+  }
+  *omax = maxv;
+  return result;
+}
+
+/**
+ * Extracts a subset of a rave data 2d field
+ * @param[in] field - the original field
+ * @param[in] ulx - the upper left x coordinate
+ * @param[in] uly - the upper left y coordinate
+ * @param[in] llx - the lower left x coordinate
+ * @param[in] lly - the lower left y coordinate
+ * @returns the subset on success otherwise NULL
+ */
+static RaveData2D_t* BeamBlockageMapInternal_subset(RaveData2D_t* field, long ulx, long uly, long llx, long lly)
+{
+  RaveData2D_t *ofield = NULL, *result = NULL;
+  long xsize = 0, ysize = 0, x = 0, y = 0;
+  RaveDataType dtype = RaveDataType_UNDEFINED;
+
+  RAVE_ASSERT((field != NULL), "field == NULL");
+
+  xsize = llx - ulx;
+  ysize = lly - uly;
+  dtype = RaveData2D_getType(field);
+
+  if (ulx < 0 || uly < 0 ||
+      llx >= RaveData2D_getXsize(field) ||
+      lly >= RaveData2D_getYsize(field) ||
+      xsize <= 0 || ysize <= 0) {
+    RAVE_WARNING0("Bad subset coorner points.");
+    goto done;
+  }
+
+  ofield = RAVE_OBJECT_NEW(&RaveData2D_TYPE);
+  if (ofield == NULL || !RaveData2D_createData(ofield, xsize, ysize, dtype)) {
+    goto done;
+  }
+
+  for (x = ulx; x < llx; x++) {
+    for (y = uly; y < lly; y++) {
+      double v = 0.0;
+      RaveData2D_getValue(field, x, y, &v);
+      RaveData2D_setValue(ofield, x-ulx, y-uly, v);
+    }
+  }
+
+  result = RAVE_OBJECT_COPY(ofield);
+done:
+  RAVE_OBJECT_RELEASE(ofield);
+  return result;
+}
+
+
+static int BeamBlockageMapInternal_getLonLatFields(PolarScan_t* scan, RaveData2D_t** lon, RaveData2D_t** lat)
+{
+  int result = 0;
+  long nrays = 0, nbins = 0;
+  long ri = 0, bi = 0;
+  double elangle = 0.0;
+  RaveData2D_t *olon = NULL, *olat = NULL;
+
+  RAVE_ASSERT((lon != NULL), "lon == NULL");
+  RAVE_ASSERT((lat != NULL), "lat == NULL");
+
+  if (scan == NULL) {
+    goto done;
+  }
+
+  nrays = PolarScan_getNrays(scan);
+  nbins = PolarScan_getNbins(scan);
+  elangle = PolarScan_getElangle(scan);
+
+  olon = RAVE_OBJECT_NEW(&RaveData2D_TYPE);
+  olat = RAVE_OBJECT_NEW(&RaveData2D_TYPE);
+
+  if (olon == NULL || olat == NULL) {
+    goto done;
+  }
+
+  if (!RaveData2D_createData(olon, nbins, nrays, RaveDataType_DOUBLE) ||
+      !RaveData2D_createData(olat, nbins, nrays, RaveDataType_DOUBLE)) {
+    goto done;
+  }
+
+  for (ri = 0; ri < nrays; ri++) {
+    for (bi = 0; bi < nbins; bi++) {
+      double lonval = 0.0, latval = 0.0;
+      if (PolarScan_getLonLatFromIndex(scan, bi, ri, &lonval, &latval)) {
+        RaveData2D_setValue(olon, bi, ri, lonval);
+        RaveData2D_setValue(olat, bi, ri, latval);
+      }
+    }
+  }
+
+  *lat = RAVE_OBJECT_COPY(olat);
+  *lon = RAVE_OBJECT_COPY(olon);
+
+  result = 1;
+done:
+  RAVE_OBJECT_RELEASE(olon);
+  RAVE_OBJECT_RELEASE(olat);
+
+  return result;
+}
 
 /**
  * Find out which maps are needed to cover given area
@@ -323,6 +489,57 @@ done:
   return result;
 }
 
+BBTopography_t* BeamBlockageMap_getTopographyForScan(BeamBlockageMap_t* self, PolarScan_t* scan)
+{
+  BBTopography_t *field = NULL, *result = NULL;
+  RaveData2D_t *lon = NULL, *lat = NULL;
+  BBTopography_t* topo = NULL;
+  double dlon_min = 0.0, dlon_max = 0.0, dlat_min = 0.0, dlat_max = 0.0;
+  double ddlon = 0.0, ddlat = 0.0;
+  long ulx = 0, uly = 0, llx = 0, lly = 0;
+
+  RAVE_ASSERT((self != NULL), "self == NULL");
+  if (scan == NULL) {
+    RAVE_ERROR0("Trying to get topography for NULL scan");
+    return NULL;
+  }
+
+  if ((topo = BeamBlockageMap_readTopography(self,
+                                             PolarScan_getLatitude(scan),
+                                             PolarScan_getLongitude(scan),
+                                             PolarScan_getMaxDistance(scan))) == NULL) {
+    goto done;
+  }
+
+  if (!BeamBlockageMapInternal_getLonLatFields(scan, &lon, &lat)) {
+    goto done;
+  }
+
+  if (!BeamBlockageMapInternal_min(lon, &dlon_min) ||
+      !BeamBlockageMapInternal_max(lon, &dlon_max) ||
+      !BeamBlockageMapInternal_min(lat, &dlat_min) ||
+      !BeamBlockageMapInternal_max(lon, &dlat_max)) {
+    goto done;
+  }
+
+  ulx = floor((dlon_min - BBTopography_getUlxmap(topo)) / BBTopography_getXDim(topo)) - 1; /* Min lon */
+  uly = floor((BBTopography_getUlymap(topo) - dlat_max) / BBTopography_getYDim(topo)) - 1; /* Max lat */
+
+  llx = ceil((dlon_max - BBTopography_getUlxmap(topo)) / BBTopography_getYDim(topo)) + 1;  /* Max Lon */
+  //lly = ceil(())
+
+  /*
+  lon_min = floor((dlon_min - BBTopography_getUlxmap(topo)) / BBTopography_getXDim(topo)) - 1;
+  lon_max = ceil((dlon_max - BBTopography_getUlxmap(topo)) / BBTopography_getXDim(topo)) + 1;
+  lat_max = floor((BBTopography_getUlymap(topo) - dlat_max) / BBTopography_getYDim(topo)) - 1;
+  lat_min = ceil((BBTopography_getUlymap(topo) - dlat_min) / BBTopography_getYDim(topo)) + 1;
+  */
+  result = RAVE_OBJECT_COPY(field);
+done:
+  RAVE_OBJECT_RELEASE(field);
+  return result;
+}
+
 /*@{ Interface functions */
 int BeamBlockageMap_setTopo30Directory(BeamBlockageMap_t* self, const char* topodirectory)
 {
index 8527e61..4bdb0fc 100644 (file)
@@ -27,6 +27,7 @@ along with beamb.  If not, see <http://www.gnu.org/licenses/>.
 #include "rave_object.h"
 #include "rave_field.h"
 #include "bbtopography.h"
+#include "polarscan.h"
 
 /**
  * Defines a beam blockage object
@@ -62,4 +63,13 @@ const char* BeamBlockageMap_getTopo30Directory(BeamBlockageMap_t* self);
  */
 BBTopography_t* BeamBlockageMap_readTopography(BeamBlockageMap_t* self, double lat, double lon, double d);
 
+/**
+ * Returns a topography that matches the scan sweep strategy. I.e. the topography
+ * for each bin/ray index.
+ * @param[in] self - self
+ * @param[in] scan - the polar scan
+ * @returns the topography with the corresponding dimensions. cols = bin and rows = rays
+ */
+BBTopography_t* BeamBlockageMap_getTopographyForScan(BeamBlockageMap_t* self, PolarScan_t* scan);
+
 #endif /* BEAMBLOCKAGEMAP_H */