diff --git a/src/stcal/outlier_detection/utils.py b/src/stcal/outlier_detection/utils.py
index 964de5fe..bc0485f3 100644
--- a/src/stcal/outlier_detection/utils.py
+++ b/src/stcal/outlier_detection/utils.py
@@ -223,7 +223,7 @@ def flag_resampled_crs(
     return mask1_smoothed & mask2
 
 
-def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
+def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio, fillval=0.0):
     """
     Resample the median data to recreate an input image based on
     the blot wcs.
@@ -236,7 +236,7 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
     median_wcs : gwcs.wcs.WCS
         The wcs for the median data.
 
-    blot_shape : list of int
+    blot_shape : tuple of int
         The target blot data shape.
 
     blot_wcs : gwcs.wcs.WCS
@@ -245,6 +245,9 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
     pix_ratio : float
         Pixel ratio.
 
+    fillval : float, optional
+        Fill value for missing data.
+
     Returns
     -------
     blotted : numpy.ndarray
@@ -259,7 +262,7 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
     log.debug("Sci shape: {}".format(blot_shape))
     log.info('Blotting {} <-- {}'.format(blot_shape, median_data.shape))
 
-    outsci = np.zeros(blot_shape, dtype=np.float32)
+    outsci = np.full(blot_shape, fillval, dtype=np.float32)
 
     # Currently tblot cannot handle nans in the pixmap, so we need to give some
     # other value.  -1 is not optimal and may have side effects.  But this is
@@ -267,7 +270,7 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
     # before a change is made.  Preferably, fix tblot in drizzle.
     pixmap[np.isnan(pixmap)] = -1
     tblot(median_data, pixmap, outsci, scale=pix_ratio, kscale=1.0,
-          interp='linear', exptime=1.0, misval=0.0, sinscl=1.0)
+          interp='linear', exptime=1.0, misval=fillval, sinscl=1.0)
 
     return outsci
 
diff --git a/tests/outlier_detection/test_utils.py b/tests/outlier_detection/test_utils.py
index c602a89c..0a59f3c2 100644
--- a/tests/outlier_detection/test_utils.py
+++ b/tests/outlier_detection/test_utils.py
@@ -118,6 +118,30 @@ def test_gwcs_blot():
     np.testing.assert_equal(blotted, median_data[:blot_shape[0], :blot_shape[1]])
 
 
+@pytest.mark.parametrize('fillval', [0.0, np.nan])
+def test_gwcs_blot_fillval(fillval):
+    # set up a very simple wcs that scales by 1x
+    output_frame = gwcs.Frame2D(name="world")
+    forward_transform = models.Scale(1) & models.Scale(1)
+
+    median_shape = (10, 10)
+    median_data = np.arange(100, dtype=np.float32).reshape((10, 10))
+    median_wcs = gwcs.WCS(forward_transform, output_frame=output_frame)
+    blot_shape = (20, 20)
+    blot_wcs = gwcs.WCS(forward_transform, output_frame=output_frame)
+    pix_ratio = 1.0
+
+    blotted = gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs,
+                        pix_ratio, fillval=fillval)
+
+    # since the blot data is larger and the wcs are equivalent the blot
+    # will contain the median data + some fill values
+    assert blotted.shape == blot_shape
+    np.testing.assert_equal(blotted[:median_shape[0], :median_shape[1]], median_data)
+    np.testing.assert_equal(blotted[median_shape[0]:, :], fillval)
+    np.testing.assert_equal(blotted[:, median_shape[1]:], fillval)
+
+
 def test_calc_gwcs_pixmap():
     # generate 2 wcses with different scales
     output_frame = gwcs.Frame2D(name="world")