Let vectors in SiliconSensor release their previous allocated memory …

…when getting a smaller target image.

include/galsim/Image.h

@@ -736,10 +736,16 @@ namespace galsim {
          *  share data with each other, but the tie to this ImageAlloc is severed.
          *
          *  This typically allocates new memory for the array.  The only
-         *  exception is if the new size is _smaller_ than currently and there
-         *  are no other views of the data.
+         *  exception is if the new size is the same as the current size and
+         *  are there are no other views of the data.  Then it just updates the
+         *  bounds to the new bounds and keeps the current array.
+         *
+         *  You can also optionally keep the current array if you are shrinking the
+         *  bounds to a smaller size, with the same limit of there not being other views.
+         *  To get this behavior, set release=false, and it will not release the allocated
+         *  memory when shrinking the size of the image.
          */
-        void resize(const Bounds<int>& new_bounds);
+        void resize(const Bounds<int>& new_bounds, bool release=true);
 
         //@{
         /**

src/Image.cpp

@@ -286,7 +286,7 @@ ImageAlloc<T>::ImageAlloc(const Bounds<int>& bounds, const T init_value) :
 }
 
 template <typename T>
-void ImageAlloc<T>::resize(const Bounds<int>& new_bounds)
+void ImageAlloc<T>::resize(const Bounds<int>& new_bounds, bool release)
 {
     if (!new_bounds.isDefined()) {
         // Then this is really a deallocation.  Clear out the existing memory.
@@ -299,9 +299,9 @@ void ImageAlloc<T>::resize(const Bounds<int>& new_bounds)
         this->_stride = 0;
         this->_ncol = 0;
         this->_nrow = 0;
-    } else if (this->_bounds.isDefined() &&
-               new_bounds.area() <= this->_nElements &&
-               this->_owner.unique()) {
+    } else if (this->_bounds.isDefined() && this->_owner.unique() &&
+               (new_bounds.area() == this->_nElements ||
+                (new_bounds.area() < this->_nElements && !release))) {
         // Then safe to keep existing memory allocation.
         // Just redefine the bounds and stride.
         this->_bounds = new_bounds;

src/Silicon.cpp

@@ -125,6 +125,8 @@ namespace galsim {
 
         _horizontalDistortions.resize(horizontalRowStride(_nx) * (_ny + 1));
         _verticalDistortions.resize(verticalColumnStride(_ny) * (_nx + 1));
+        _horizontalDistortions.shrink_to_fit();
+        _verticalDistortions.shrink_to_fit();
 
         for (int index=0; index < nv1*_nx*_ny; index++) {
             int n1 = index % nv1;
@@ -1007,6 +1009,8 @@ namespace galsim {
     {
         _horizontalBoundaryPoints.resize(horizontalRowStride(nx) * (ny+1));
         _verticalBoundaryPoints.resize(verticalColumnStride(ny) * (nx+1));
+        _horizontalBoundaryPoints.shrink_to_fit();
+        _verticalBoundaryPoints.shrink_to_fit();
 
         // fill in horizontal boundary points from emptypoly
         int i = 0;
@@ -1037,6 +1041,8 @@ namespace galsim {
 
         _pixelInnerBounds.resize(nx * ny);
         _pixelOuterBounds.resize(nx * ny);
+        _pixelInnerBounds.shrink_to_fit();
+        _pixelOuterBounds.shrink_to_fit();
         for (int k = 0; k < (nx * ny); k++) {
             updatePixelBounds(nx, ny, k, _pixelInnerBounds.data(),
                               _pixelOuterBounds.data(),