Repository: spacetelescope/jwst

Issue:
Refactor model_blender to take in only necessary metadata
_Issue [JP-3701](https://jira.stsci.edu/browse/JP-3701) was created on JIRA by [Ned Molter](https://jira.stsci.edu/secure/ViewProfile.jspa?name=emolter):_

The model_blender suite, which is called during calwebb_image3 by default (inside resample.resample_many_to_one), expects to take in an in-memory list of datamodels as input.  This poses a problem for the promised memory savings migration to ModelLibrary; we want to avoid loading multiple models into memory at once.  There are two options for how this could be fixed in the long-term:
 * Make model_blender expect a ModelLibrary, and load models into memory as needed within that step using ModelLibrary's architecture that is designed to handle this.  Looking into the source code a little bit, it seems that the list of input models is accessed in several different places, which means that a borrow...shelve (load...save) operation would be needed several times.  This would be fairly straightforward to do, but might be slow
 * Make model_blender take a list of metadata objects as input.  It appears that the data are not needed anywhere in the step, and just a list of metadata objects would obviously take up much less memory.  In that case only one borrow...shelve (load...save) operation would be needed for each model to append its model.meta to a list.  This solution takes more refactoring of model_blender but in my opinion would be the preferred way to go.

Setting this to low priority because there is a workaround (suggested by [Melanie Clarke](https://jira.stsci.edu/secure/ViewProfile.jspa?name=mclarke)), which is to set all the data arrays of each model to be empty before passing them as an in-memory list to model_blender.  This workaround is not yet tested but seems likely to work without bad consequences in the short term.

 

Discussion:
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) Report
Attention: Patch coverage is `92.56757%` with `11 lines` in your changes missing coverage. Please review.
> Project coverage is 60.73%. Comparing base [(`8381a26`)](https://app.codecov.io/gh/spacetelescope/jwst/commit/8381a261d0b9a51d0ce2e7ac693032c76efac697?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) to head [(`5e32fb5`)](https://app.codecov.io/gh/spacetelescope/jwst/commit/5e32fb58a60ca7b0c5146ecd75a1909aecef171d?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope).
> Report is 3 commits behind head on master.

| [Files with missing lines](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?dropdown=coverage&src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) | Patch % | Lines |
|---|---|---|
| [jwst/model\_blender/rules.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fmodel_blender%2Frules.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9tb2RlbF9ibGVuZGVyL3J1bGVzLnB5) | 87.87% | [4 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/model\_blender/blender.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fmodel_blender%2Fblender.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9tb2RlbF9ibGVuZGVyL2JsZW5kZXIucHk=) | 94.59% | [2 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/cube\_build/ifu\_cube.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fcube_build%2Fifu_cube.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9jdWJlX2J1aWxkL2lmdV9jdWJlLnB5) | 0.00% | [1 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/model\_blender/schemautil.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fmodel_blender%2Fschemautil.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9tb2RlbF9ibGVuZGVyL3NjaGVtYXV0aWwucHk=) | 96.00% | [1 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/model\_blender/tablebuilder.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fmodel_blender%2Ftablebuilder.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9tb2RlbF9ibGVuZGVyL3RhYmxlYnVpbGRlci5weQ==) | 97.36% | [1 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/pipeline/calwebb\_coron3.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fpipeline%2Fcalwebb_coron3.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9waXBlbGluZS9jYWx3ZWJiX2Nvcm9uMy5weQ==) | 0.00% | [1 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |
| [jwst/resample/resample.py](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&filepath=jwst%2Fresample%2Fresample.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope#diff-andzdC9yZXNhbXBsZS9yZXNhbXBsZS5weQ==) | 85.71% | [1 Missing :warning: ](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope) |

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff             @@
##           master    #8759      +/-   ##
==========================================
+ Coverage   60.63%   60.73%   +0.09%     
==========================================
  Files         372      374       +2     
  Lines       38372    38262     -110     
==========================================
- Hits        23267    23238      -29     
+ Misses      15105    15024      -81     
```

</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/spacetelescope/jwst/pull/8759?dropdown=coverage&src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=spacetelescope).

braingram: I fixed the 2 bugs found while working on this PR (the ones mentioned in the PR description).

This results in many regtest "failures":
https://plwishmaster.stsci.edu:8081/blue/organizations/jenkins/RT%2FJWST-Developers-Pull-Requests/detail/JWST-Developers-Pull-Requests/1698/tests

with "[stable-deps] test_miri_lrs_slit_spec3[default_wcs-x1d] – jwst.regtest.test_miri_lrs_slit_spec3" being a good example showing both bugs. With this PR the fitsdiff fails with:
```
     Headers contain differences:
       Headers have different number of cards:
        a: 256
        b: 257
       Extra keyword 'OPMODE' in b: 'NONE'
       Keyword TIME-OBS has different values:
          a> 23:07:34.995
           ?            ^
          b> 23:07:34.994999
           ?   
```
The `TIME-OBS` in the result (a) is more accurate than truth (b) and model blender no longer creates the default `OPMODE` (which doesn't exist in the input, blended models).

There are also many unrelated regtest failures in this run due removal of `R_DRZPAR`, what looks like a crds update and existing failures.
braingram: Thanks for the suggestions and review @melanieclarke I think I addressed all your comments but please let me know if there are any I missed and if you have any follow-up comments from the changes.
melanieclarke: > Are the regtests in a state where the changes introduced by this PR can be okified or are there current failures?

They should be pretty clear today, when you're ready to re-run.

Patch:
diff --git a/CHANGES.rst b/CHANGES.rst
index a17ab3914..4108daff1 100644
--- a/CHANGES.rst
+++ b/CHANGES.rst
@@ -95,6 +95,11 @@ master_background
 - Either of ``"background"`` or ``"bkg"`` in slit name now defines the slit
   as a background slit, instead of ``"bkg"`` only. [#8600]
 
+model_blender
+-------------
+
+- Allow incremental blending of models. [#8759]
+
 mrs_imatch
 ----------
 
diff --git a/docs/jwst/model_blender/blender_rules.rst b/docs/jwst/model_blender/blender_rules.rst
index edafb7061..55d5e6d67 100644
--- a/docs/jwst/model_blender/blender_rules.rst
+++ b/docs/jwst/model_blender/blender_rules.rst
@@ -5,39 +5,15 @@ Model Blender Rules
 
 Blending models relies on rules to define how to evaluate all the input values
 for a model attribute in order to determine the final output value. These rules
-then get specified in the model schema for each attribute.
+are derived from the model schema for each attribute.
 
-The rules get interpreted and applied as list or array operations that work on
-the set of input values for each attribute.  The full set of pre-defined rules
-includes
+The rules are applied to a collection of all values to be blended (for example
+if blending models that have different exposure times, the blended exposure
+time will be the sum of all exposure times).
 
-.. code-block:: python
+Supported rule names and corresponding functions are listed in
+`jwst.model_blender.rules.RULE_FUNCTIONS`.
 
-    import numpy as np
-    # translation dictionary for function entries from rules files
-    blender_funcs = {'first': first,
-                      'last': last,
-                      'float_one': float_one,
-                      'int_one': int_one,
-                      'zero': zero,
-                      'multi': multi,
-                      'multi?': multi1,
-                      'mean': np.mean,
-                      'sum': np.sum,
-                      'max': np.max,
-                      'min': np.min,
-                      'stddev': np.std,
-                      'mintime': mintime,
-                      'maxtime': maxtime,
-                      'mindate': mindate,
-                      'maxdate': maxdate,
-		      'mindatetime': mindatetime,
-		      'maxdatetime': maxdatetime}
-
-The rules that should be referenced in the model schema definition are the
-keys defined for `jwst.model_blender.blender_rules.blender_funcs` listed
-above.  This definition illustrates how several rules are simply interfaces for
-numpy array operations, while others are defined internally to `model_blender`.
-
-.. automodapi:: jwst.model_blender.blendrules
+.. automodapi:: jwst.model_blender.rules
   :no-inheritance-diagram:
+  :include-all-objects:
diff --git a/docs/jwst/model_blender/index.rst b/docs/jwst/model_blender/index.rst
index fd76d52ea..ae9022214 100644
--- a/docs/jwst/model_blender/index.rst
+++ b/docs/jwst/model_blender/index.rst
@@ -10,4 +10,3 @@ Model Blender
    blender_rules.rst
 
 .. automodapi:: jwst.model_blender
-   :skip: "OrderedDict"
diff --git a/docs/jwst/model_blender/main.rst b/docs/jwst/model_blender/main.rst
index fc970dba3..97b8bb561 100644
--- a/docs/jwst/model_blender/main.rst
+++ b/docs/jwst/model_blender/main.rst
@@ -3,6 +3,14 @@
 Role of Model Blender
 ======================
 
+Model Blender performs 2 operations:
+
+  - Accumulate, blend and set metadata when combining several models (
+    as is done during resampling) using combination rules defined in the datamodel
+    schemas.
+  - Record the metadata of input models in a table added to the
+    combined model.
+
 One problem with combining data from multiple exposures stems from not being able
 to keep track of what kind of data was used to create the final product.  The
 final product only reports one value for each of the metadata attributes from the
@@ -10,7 +18,7 @@ model schema used to describe the science data, where each of multiple inputs ma
 have had different values for each attribute.  The model_blender package solves
 this problem by allowing the user to define rules that can be used to determine a
 final single value from all the input values for each attribute, using separate
-rules for each attribute as appropriate.  This package also creates a FITS binary
+rules for each attribute as appropriate.  This package also creates a
 table that records the input attribute values for all the input models used to
 create the final product, allowing the user to select what attributes to keep in
 this table.
@@ -27,86 +35,98 @@ This code works by
 
 Using model_blender
 ===================
-The model blender package requires
-
-  - all the input models be available
-  - the output product has already been generated
-
-Both the input models and output product could be provided as either a datamodel
-instance in memory or as the name of a FITS file on disk.  The primary advantage
-to working strictly in-memory with datamodel instances comes from minimizing the
-amount of disk I/O needed for this operation which can result in significantly
-more efficient (read that: faster) processing.
 
-.. note::
 
-  The generated output model will be considered to contain a default
-  (or perhaps even empty) set of :ref:`metadata` based on some
-  model defined in :ref:`stdatamodels<stdatamodels:data-models>`
-  This metadata will be replaced
-  **in-place** when running :ref:`blender_api`.
+Using blendmodels
+-----------------
 
 The simplest way to run model blender only requires calling a single interface:
 
 .. code-block:: python
 
   from jwst.model_blender import blendmeta
-  blendmeta.blendmodels(product, inputs=input_list)
+  blendmeta.blendmodels(product, input_list)
 
 where
 
-  - `product`: the datamodel (or FITS filename) for the already combined product
-  - `input_list`: list of input datamodels or FITS filenames for all inputs used
-    to create the `product`
+  - `product`: the datamodel for the already combined product
+  - `input_list`: list of input datamodels used to create the `product`
 
 
-The output product will end up with new metadata attribute values and a new HDRTAB
-FITS binary table extension in the FITS file when the product gets saved to disk.
+The output product will end up with new metadata attribute values and a new ``hdrtab``
+attribute that will produce a ``HDRTAB`` table extension when the product is saved
+to disk.
+
+Using ModelBlender
+------------------
+
+Alternatively, metadata blending can be done incrementally by using `ModelBlender`
+to accumulate metadata from several models and then blending the result into an
+output model.
+
+.. code-block:: python
+
+   from jwst.model_blender import ModelBlender
+   blender = ModelBlender()
+   for model in input_models:
+       blender.accumulate(model)
+   blender.finalize_model(product)
+
+This produces a ``product`` identical to a call to `blendmodels` described above.
 
 
 Customizing the behavior
 ========================
-By default, `blendmodels` will not write out the updated product model to disk.
-This allows the user or calling program to revise or apply data-specific logic
-to redefine the output value for any of the output product's metadata attributes.
-For example, when combining multiple images, the WCS information does not represent
-any combination of the input WCS attributes.  Instead, the user can have
-their own processing code replace the *blended* WCS attributes with one that was
-computed separately using a complex, accurate algorithm.  This is, in fact, what
-the resample step does to create the final resampled output product whenever it is
-called by steps in the JWST pipeline.
-
-Additional control over the behavior of model_blender comes from editing the
-schema for the input datamodels where the rules for each attribute are defined.
-A sample definition from the core schema demonstrates the basic syntax used for
-any model blending definitions::
-
-          time_end:
-            title: UTC time at end of exposure
-            type: string
-            fits_keyword: TIME-END
-            blend_rule: last
+
+
+Combined model metadata
+-----------------------
+
+The metadata assigned to the combined (output) model will be determined from:
+
+  - the first combined model metadata
+  - any metadata with a ``blend_rule``
+
+Looking at this sample from the core schema::
+
+          start_time:
+            title: "[d] exposure start time in MJD"
+            type: number
+            fits_keyword: EXPSTART
             blend_table: True
+            blend_rule: min
+
+Since ``start_time`` uses ``blend_rule: min`` the ``start_time`` for the combined
+model will be the minimum of all input models. If ``blend_rule`` was not defined, the
+``start_time`` for the combined model would be the ``start_time`` of the first model.
+
+.. warning::
+
+   All blended models must be of the same type.
+
+.. warning::
+
+   Since model blender copies the metadata of the first model, all metadata
+   from that model will be added to combined model. This is true even for
+   metadata that is not defined in the schema.
+
+.. warning::
+
+   Metadata nested in arrays will not be blended. For example,
+   ``meta.background.polynomial_info`` contains a list/array of
+   dictionaries. Model blender will not blend this attribute.
+
+Input model metadata Table
+--------------------------
+
+As seen in the above example, ``start_time`` defines ``blend_table: True``.
+This tells model blender to include this attribute in the table added to the
+combined model that lists attributes from input models. If the attribute
+has a ``fits_keyword``, the resulting column will use the keyword for the
+column name (if not defined, the attribute name will be used). So for
+the above example, the resulting table will contain an ``EXPSTART`` column
+(and not a ``meta.exposure.start_time`` column) with a row listing
+the ``start_time`` for each input model.
 
-Any attribute without an entry for `blend_rule` will use the default rule of
-`first` which selects the first value from all inputs in the order provided as the
-final output value.  Any attribute with a `blend_table` rule will insure that
-the specific attribute will be included in the output HDRTAB binary table appended
-to the product model when it gets written out to disk as a FITS file.
-
-The full set of rules included in the package are described in
-:ref:`blender_rules` and include common list/array operations such as
-(but not limited to):
-
-  - minimum
-  - maximum
-  - first
-  - last
-  - mean
-  - zero
-
-These can then be used to customize the output value for any given attribute
-should the rule provided by default with the schema installed with the
-JWST environment not be correct for the user's input data.  The user can simply
-edit the schema definition installed in their JWST environment to apply custom
-rules for blending the data being processed.
+If an input model is missing an attribute a ``nan`` will be stored in the
+corresponding table cell.
diff --git a/jwst/cube_build/ifu_cube.py b/jwst/cube_build/ifu_cube.py
index 3cd627338..ed034949d 100644
--- a/jwst/cube_build/ifu_cube.py
+++ b/jwst/cube_build/ifu_cube.py
@@ -2362,10 +2362,14 @@ class IFUCubeData():
     # ********************************************************************************
     def blend_output_metadata(self, IFUCube):
         """Create new output metadata based on blending all input metadata."""
-        # Run fitsblender on output product
-        output_file = IFUCube.meta.filename
-        blendmeta.blendmodels(IFUCube, inputs=self.input_models_this_cube,
-                              output=output_file)
+        # Run blendmodels on output product
+        blendmeta.blendmodels(
+            IFUCube,
+            self.input_models_this_cube,
+            ignore=[
+                "meta.filename",
+            ],
+        )
 
 
 class IncorrectInput(Exception):
diff --git a/jwst/model_blender/__init__.py b/jwst/model_blender/__init__.py
index e69de29bb..a2128990d 100644
--- a/jwst/model_blender/__init__.py
+++ b/jwst/model_blender/__init__.py
@@ -0,0 +1,5 @@
+from .blender import ModelBlender
+from .blendmeta import blendmodels
+
+
+__all__ = ['ModelBlender', 'blendmodels']
diff --git a/jwst/model_blender/_schemautil.py b/jwst/model_blender/_schemautil.py
new file mode 100644
index 000000000..be84a9c4e
--- /dev/null
+++ b/jwst/model_blender/_schemautil.py
@@ -0,0 +1,67 @@
+from stdatamodels import schema as dm_schema
+
+def parse_schema(schema):
+    """
+    Parse an ASDF schema for model blending instructions.
+
+    Parameters
+    ----------
+    schema : dict
+        Dictionary containing an ASDF schema.
+
+    Returns
+    -------
+    attr_to_column: dict
+        one-to-one mapping of metadata attributes
+        (full dotted-paths) to table columns. For example
+        {'meta.observation.time': 'TIME-OBS'} for storing
+        the 'meta.observation.time' attributes in a 'TIME-OBS'
+        column.
+
+    attr_to_rule: dict
+        mapping of metadata attribute to blend rule
+        type. For example {'meta.observation.time': 'mintime'}
+        will combine all 'meta.observation.time' attributes
+        using the 'mintime' rule.
+
+    schema_ignores: list
+        list of attributes that will not be blended
+    """
+    def callback(subschema, path, combiner, ctx, recurse):
+        if len(path) <= 1:
+            return  # ignore top-level (non-meta) attributes
+        if path[0] != 'meta':
+            return  # ignore non-meta attributes
+        if 'items' in path:
+            return  # ignore attributes in arrays
+        if subschema.get('properties'):
+            return  # ignore ObjectNodes
+
+        # strip trailing path if there's a combiner
+        for combiner in ['anyOf', 'oneOf']:
+            if combiner in path:
+                path = path[:path.index(combiner)]
+                break
+
+        # construct the metadata attribute path
+        attr = '.'.join(path)
+
+        if subschema.get('type') == 'array':
+            ctx['ignores'].append(attr)
+            return  # ignore ListNodes
+
+        # if 'blend_rule' is defined, make a 'blend'
+        if 'blend_rule' in subschema:
+            ctx['blends'][attr] = subschema['blend_rule']
+
+        # if 'blend_table' is defined (and truthy), add a column
+        if subschema.get('blend_table'):
+            ctx['columns'][attr] = subschema.get('fits_keyword', attr)
+
+    ctx = {
+        'columns': {},
+        'blends': {},
+        'ignores': [],
+    }
+    dm_schema.walk_schema(schema, callback, ctx)
+    return ctx['columns'], ctx['blends'], ctx['ignores']
diff --git a/jwst/model_blender/_tablebuilder.py b/jwst/model_blender/_tablebuilder.py
new file mode 100644
index 000000000..e677d0699
--- /dev/null
+++ b/jwst/model_blender/_tablebuilder.py
@@ -0,0 +1,118 @@
+import numpy as np
+
+
+def _convert_dtype(value):
+    """Convert numarray column dtype into YAML-compatible format description"""
+    if 'U' in value:
+        # working with a string description
+        str_len = int(value[value.find('U') + 1:])
+        new_dtype = ['ascii', str_len]
+    elif value == 'bool':
+        # cast all bool to int8 to avoid issues on write
+        new_dtype = 'int8'
+    else:
+        new_dtype = str(value)
+
+    return new_dtype
+
+
+def table_to_schema(table):
+    """
+    Construct a schema for a table.
+
+    Convert a "table" (a structured ndarray) to a stdatamodels
+    sub-schema that will allow the "table" to be stored to a fits
+    extension HDRTAB.
+
+    Parameters
+    ----------
+    table: ndarray
+        The structured array containing the data (and datatype).
+
+    Returns
+    -------
+    subschema: dict
+        stdatamodels for the "table" datatype.
+    """
+    return {
+        'title': 'Combined header table',
+        'fits_hdu': 'HDRTAB',
+        'datatype': [
+            {
+                'name': col_name,
+                'datatype': _convert_dtype(str(table.dtype[col_name])),
+            }
+            for col_name in table.dtype.fields
+        ],
+    }
+
+
+class TableBuilder:
+    """
+    Class to build a metadata table.
+
+    Used to incrementally build a metadata "table" (a numpy
+    structured array) containing metadata from several models.
+
+    >>> tb = TableBuilder({"meta.filename": "FN"})
+    >>> tb.header_to_row({"meta.filename": "foo.fits"})
+    >>> tb.build_table()
+    rec.array([('foo.fits',)],
+         dtype=[('FN', '<U8')])
+    """
+    def __init__(self, attr_to_column):
+        """
+        Create a new `TableBuilder`.
+
+        Parameters
+        ----------
+        attr_to_column: dict
+            A one-to-one mapping of attribute names (as
+            dotted paths like "meta.filename") to column
+            names (strings).
+        """
+        self.attr_to_column = attr_to_column
+        self.columns = {col: [] for col in self.attr_to_column.values()}
+        if len(attr_to_column) != len(self.columns):
+            raise ValueError(f"Invalid attr_to_column, mapping is not 1-to-1: {attr_to_column}")
+
+    def header_to_row(self, header):
+        """
+        Add metadata in header to the table.
+
+        This function will add a complete row for each
+        header. If header is missing a required attribute
+        (as defined in the mapping provided to
+        `TableBuilder.__init__`) the column with this
+        missing value for this row will contain a ``nan``.
+
+        Parameters
+        ----------
+        header: dict
+            Often produced from ``Datamodel.to_flat_dict``.
+        """
+        row = {}
+        for attr in self.attr_to_column:
+            if attr in header:
+                row[attr] = header[attr]
+        self._add_row(row)
+
+    def _add_row(self, row):
+        for attr, col in self.attr_to_column.items():
+            self.columns[col].append(row[attr] if attr in row else np.nan)
+
+    def build_table(self):
+        """
+        Build the final "table".
+
+        Returns
+        -------
+        table: numpy.ndarray
+            Structured array containing fields with datatypes
+        """
+        arrays = []
+        table_dtype = []
+        for col, items in self.columns.items():
+            arrays.append(np.array(items))
+            table_dtype.append((col, arrays[-1].dtype))
+        return np.rec.fromarrays(arrays, dtype=table_dtype)
diff --git a/jwst/model_blender/blender.py b/jwst/model_blender/blender.py
index bb4afe487..732eb6455 100644
--- a/jwst/model_blender/blender.py
+++ b/jwst/model_blender/blender.py
@@ -1,254 +1,158 @@
-# Copyright (C) 2010 Association of Universities for Research in Astronomy(AURA)
-#
-# Redistribution and use in source and binary forms, with or without
-# modification, are permitted provided that the following conditions are met:
-#
-#     1. Redistributions of source code must retain the above copyright
-#       notice, this list of conditions and the following disclaimer.
-#
-#     2. Redistributions in binary form must reproduce the above
-#       copyright notice, this list of conditions and the following
-#       disclaimer in the documentation and/or other materials provided
-#       with the distribution.
-#
-#     3. The name of AURA and its representatives may not be used to
-#       endorse or promote products derived from this software without
-#       specific prior written permission.
-#
-# THIS SOFTWARE IS PROVIDED BY AURA ``AS IS'' AND ANY EXPRESS OR IMPLIED
-# WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
-# MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
-# DISCLAIMED. IN NO EVENT SHALL AURA BE LIABLE FOR ANY DIRECT, INDIRECT,
-# INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
-# BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS
-# OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
-# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR
-# TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE
-# USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
-# DAMAGE.
-
-import numpy as np
-from numpy import ma
-
-from stdatamodels.jwst.datamodels import JwstDataModel
-from stdatamodels.jwst import datamodels
-
-
-class _KeywordMapping:
-    """
-    A simple class to verify and store information about each mapping
-    entry.
-    """
-
-    def __init__(self, src_kwd, dst_name, agg_func=None, error_type="ignore",
-                 error_value=np.nan):
-        if not isinstance(src_kwd, str):
-            raise TypeError(
-                "The source keyword name must be a string")
-
-        if not isinstance(dst_name, str):
-            raise TypeError(
-                "The destination name must be a string")
-
-        if agg_func is not None:
-            try:
-                for i in agg_func:
-                    if not hasattr(i, '__call__'):
-                        raise TypeError(
-                            "The aggregating function must be a callable " +
-                            "object, None or a sequence of callables")
-                self.agg_func_is_sequence = True
-            except TypeError:
-                if not hasattr(agg_func, '__call__'):
-                    raise TypeError(
-                        "The aggregating function must be a callable object, "
-                        "None or a sequence of callables")
-                self.agg_func_is_sequence = False
-
-        if error_type not in ('ignore', 'raise', 'constant'):
-            raise ValueError(
-                "The error type must be either 'ignore', 'raise' or 'constant'")
-
-        self.src_kwd = src_kwd
-        self.dst_name = dst_name
-        self.agg_func = agg_func
-        self.error_type = error_type
-        self.error_value = error_value
-
-
-def metablender(input_models, spec):
-    """
-    Given a list of datamodels, aggregate metadata attribute values and
-    create a table made up of values from a number of metadata instances,
-    according to the given specification.
-
-    **Parameters:**
-
-    - *input_models* is a sequence where each element is either:
-
-      - a `datamodels.JwstDataModel` instance or sub-class
-
-      - a string giving the *filename* for the input_model
-
-    - *spec* is a list defining which keyword arguments are to be
-      aggregated and how.  Each element in the list should be a
-      sequence with 2 to 5 elements of the form:
-
-        (*src_keyword*, *dst_name*, *function*, *error_type*, *error_value*)
-
-      - *src_keyword* is the keyword to pull values from.  It is
-        case-insensitive.
+from astropy.io import fits
+from stdatamodels import fits_support
 
-      - *dst_name* is the name to use as a dictionary key or column
-        name for the destination values.
+from .rules import make_blender
+from ._schemautil import parse_schema
+from ._tablebuilder import TableBuilder, table_to_schema
 
-      - *function* (optional).  If function is not None, the values
-        from the source are aggregated and returned in the
-        *aggregate_dict*.  If function is None (or the tuple contains
-        only 2 elements), all values are stored as a column with the
-        name *dst_name* in the result *table*.
 
-        If not None, *function* should be a callable object that takes
-        a sequence of values and returns an aggregate result.  If the
-        function returns None, no values will be added to the
-        aggregate dictionary.  There are many functions in Numpy that
-        are directly useful as an aggregating function, for example:
+__all__ = ['ModelBlender']
 
-          - mean: `numpy.mean`
 
-          - median: `numpy.median`
-
-          - maximum: `numpy.max`
-
-          - minimum: `numpy.min`
-
-          - sum: `numpy.sum`
-
-          - standard deviation: `numpy.std`
-
-        Lambda functions are also often useful:
-
-          - first: ``lambda x: x[0]``
-
-          - last: ``lambda x: x[-1]``
-
-        Additionally, *function* may be a tuple, where each member is
-        itself a callable object.  The result will be a tuple
-        containing results from each of the given functions.  For
-        instance, to aggregate a range of values, i.e. both the
-        minimum and maximum values, use the following as *function*:
-        ``(numpy.min, numpy.max)``.
-
-      - *error_type* (optional) defines how missing or syntax-errored
-        values are handled.  It may be one of the following:
-
-        - 'ignore': missing or unparsable values are ignored.  They
-          are not included in the list of values passed to the
-          aggregating function.  In the result *table*, missing values
-          are masked out.
-
-        - 'raise': missing or unparsable values raise a `ValueError`
-          exception.
+class ModelBlender:
+    """
+    Class to "blend" metadata from several datamodels.
 
-        - 'constant': missing or unparsable values are replaced with a
-          constant, given by the *error_value* field.
+    Tht output "blended" model will contain:
 
-      - *error_value* (optional) is the constant value to be used for
-        missing or unparsable values when *error_type* is set to
-        'constant'.  When not provided, it defaults to `NaN`.
+        - metadata for a combined model
+        - a table with metadata of each datamodel
 
-    **Returns:**
+    Input models can be added to the blender using `ModelBlender.accumulate`
+    and the output/combined model updated using `ModelBlender.finalize_model`.
 
-    A 2-tuple of the form (*aggregate_dict*, *table*) where:
+    All input/accumulated models must be of the same type.
 
-    - *aggregate_dict* is a dictionary of where the keys come from
-      *dst_name* and the values are the aggregated values as run_KeywordMapping
-      through *function*.
+    >>> blender = ModelBlender()
+    >>> blender.accumulate(input_model_a)  # doctest: +SKIP
+    >>> blender.accumulate(input_model_b)  # doctest: +SKIP
+    >>> blender.finalize_model(combined_model)  # doctest: +SKIP
 
-    - *table* is a masked Numpy structured array where the column
-      names come from *dst_name* and the column contains the values
-      from *src_keyword* for all of the given headers.  Missing values
-      are masked out.
     """
-    mappings = [_KeywordMapping(*x) for x in spec]
-    data = [[] for x in spec]
-    data_masks = [[] for x in spec]
-
-    # Read in data
-    for model in input_models:
-        if not isinstance(model, JwstDataModel):
-            if not isinstance(model, str):
-                raise TypeError(
-                    "Each entry in the headers list must be either a " +
-                    "datamodels.JwstDataModel instance or a filename (str)")
-            model = datamodels.open(model)
-        header = model.to_flat_dict()
-        filename = header['meta.filename']
-        for i, mapping in enumerate(mappings):
-            if mapping.src_kwd in header:
-                value = model[mapping.src_kwd]
-            elif mapping.error_type == 'raise':
+    def __init__(self, blend_ignore_attrs=None):
+        """
+        Create a new `ModelBlender`.
+
+        Parameters
+        ----------
+        blend_ignore_attrs: list or None
+            A list of metadata attributes to ignore during blending.
+            These attributes will not be set on the output/combined.
+            These attributes must be strings containing the dotted
+            path of each attribute (for example "meta.filename").
+            (Note that "meta.wcs" will always be ignored).
+        """
+        self._model_type = None
+        self._first_header_meta = None
+        self._blenders = None
+        self._table_builder = None
+        self._blend_ignore_attrs = ['meta.wcs']
+        if blend_ignore_attrs is not None:
+            self._blend_ignore_attrs.extend(blend_ignore_attrs)
+
+    def accumulate(self, model):
+        """
+        Blend metadata for model.
+
+        Process model adding its metadata to the blended
+        metadata and the metadata table.
+
+        Parameters
+        ----------
+        model: `jwst.datamodels.JwstDataModel`
+            The datamodel to blend.
+        """
+        if self._first_header_meta is None:
+            self._model_type = type(model)
+            # search the schema for other metadata to "blend" and to add to the table
+            attr_to_columns, attr_to_blend_rules, schema_ignores  = parse_schema(model.schema)
+
+            # update ignores list for items in schema that can't be blended
+            self._blend_ignore_attrs.extend(schema_ignores)
+
+            # capture the entire contents of the first model metadata
+            self._first_header_meta = {}
+            for attr, v in model.to_flat_dict(include_arrays=False).items():
+                if not attr.startswith('meta'):
+                    continue
+                if any((attr.startswith(i) for i in self._blend_ignore_attrs)):
+                    continue
+                self._first_header_meta[attr] = v
+
+            # make "blenders" for the metadata with special rules
+            self._blenders = {}
+            for attr, rule in attr_to_blend_rules.items():
+                if rule == 'first':
+                    continue
+                if any((attr.startswith(i) for i in self._blend_ignore_attrs)):
+                    continue
+                self._blenders[attr] = make_blender(rule)
+
+            # make a table builder using the mapping from the schema
+            self._table_builder = TableBuilder(attr_to_columns)
+        else:
+            if type(model) != self._model_type:  # noqa: E721
                 raise ValueError(
-                    "%s is missing keyword '%s'" %
-                    (filename, mapping.src_kwd))
-            elif mapping.error_type == 'constant':
-                value = mapping.error_value
-            else:
-                value = None
-
-            if mapping.agg_func is None:
-                if value is None:
-                    data[i].append(np.nan)
-                    data_masks[i].append(True)
-                else:
-                    data[i].append(value)
-                    data_masks[i].append(False)
-            else:
-                if value is not None:
-                    data[i].append(value)
-
-    # Aggregate data into dictionary
-    results = {}
-    for i, mapping in enumerate(mappings):
-        if data[i] == []:
-            result = None
-            continue
-        if mapping.agg_func is not None:
-            if mapping.agg_func_is_sequence:
-                result = []
-                for func in mapping.agg_func:
-                    result.append(func(data[i]))
-                result = tuple(result)
-            else:
-                result = mapping.agg_func(data[i])
-            if result is not None:
-                results[mapping.dst_name] = result
-
-    # Aggregate data into table
-    dtype = []
-    arrays = []
-
-    # Use Numpy to "guess" a data type for each of the columns
-    for i, mapping in enumerate(mappings):
-        if mapping.agg_func is None:
-            array = np.array(data[i])
-            if np.issubdtype(np.int32, array.dtype):
-                # see about recasting as int32
-                if not np.any(array / (2**31 - 1) > 1.):
-                    array = array.astype(np.int32)
-            dtype.append((mapping.dst_name, array.dtype))
-            arrays.append(array)
-
-    if len(dtype):
-        # Combine the columns into a structured array
-        table = ma.empty((len(input_models),), dtype=dtype)
-        j = 0
-        for i, mapping in enumerate(mappings):
-            if mapping.agg_func is None:
-                table.data[mapping.dst_name] = arrays[j]
-                table.mask[mapping.dst_name] = data_masks[i]
-                j += 1
-    else:
-        table = np.empty((0,))
-
-    return results, table
+                    f"model of type {type(model)} does not match previous type({self._model_type}). "
+                    "ModelBlender only supports blending models of the same model type."
+                )
+
+        # convert the model to a flat header
+        header = model.to_flat_dict(include_arrays=False)
+
+        # add the header to the table
+        self._table_builder.header_to_row(header)
+
+        # and perform any special blending
+        for attr, blender in self._blenders.items():
+            if attr in header:
+                blender.accumulate(header[attr])
+
+    def _finalize_metadata(self):
+        # start with the entire contents of the first model
+        meta = self._first_header_meta.copy()
+        for attr, blender in self._blenders.items():
+            meta[attr] = blender.finalize()
+        return meta
+
+    def _finalize_table(self):
+        return self._table_builder.build_table()
+
+    def finalize_model(self, model):
+        """
+        Update model with the blend results.
+
+        Add blended metadata and the accumulated metadata table to
+        the provided datamodel. The update process involves:
+
+            - setting the model metadata to the blended metadata values
+            - adding an "hdrtab" attribute (containing the metadata table)
+            - updating the model schema to save "hdrtab"
+
+        The provided model will be updated in-place.
+
+        Parameters
+        ----------
+        model: `jwst.datamodels.JwstDataModel`
+            A datamodel that will have its metadata set
+            to the blended metadata and have the metadata
+            table assigned to the "hdrtab" attribute.
+        """
+        # update metadata of the output model based on the results
+        # of the "blenders"
+        for attr, val in self._finalize_metadata().items():
+            try:
+                model[attr] = val
+            except KeyError:
+                # Ignore keys that are in the asdf tree but not in the schema
+                pass
+
+        # patch the table into the output model and the schema
+        table = self._finalize_table()
+        schema = table_to_schema(table)
+        model.add_schema_entry('hdrtab', schema)
+
+        # because astropy will silently mangle boolean columns on write
+        # we roundtrip the data through a BinTableHDU here to allow
+        # stdatamodels to correct the data in a way that won't result in mangling
+        model.hdrtab = fits_support.from_fits_hdu(fits.BinTableHDU.from_columns(table), schema)
diff --git a/jwst/model_blender/blendmeta.py b/jwst/model_blender/blendmeta.py
index 9d6888614..d9a5434e0 100644
--- a/jwst/model_blender/blendmeta.py
+++ b/jwst/model_blender/blendmeta.py
@@ -3,328 +3,36 @@
     This module will create a new metadata instance and table from a
     list of input datamodels or filenames.
 """
-from collections import OrderedDict
 
-from astropy.io import fits
+from .blender import ModelBlender
 
-from stdatamodels import fits_support
-from stdatamodels import schema as dm_schema
-from stdatamodels.jwst import datamodels
 
-from jwst.datamodels import ModelContainer
-
-from .. import associations
-
-from .blendrules import KeywordRules
-
-
-EMPTY_LIST = [None, '', ' ', 'INDEF', 'None']
-
-__doctest_skip__ = ['blendmodels']
+__all__ = ['blendmodels']
 
 # Primary functional interface for the code
 
-
-def blendmodels(product, inputs=None, output=None, ignore=None, verbose=False):
-    """
-    Run main interface for blending metadata from multiple models.
-
-    Blend models that went into creating the original drzfile into a
-    new metadata instance with a table that contains attribute values from
-    all input datamodels.
-
-    The product will be used to determine the names of the input models,
-    should no filenames be provided in the 'inputs' parameter.
-
-    The product will be updated 'in-place' with the new metadata attributes
-    and FITS BinTableHDU table.  The blended FITS table, with extname=HDRTAB,
-    has 1 column for each metadata attribute recorded from the input models,
-    one row for each input model, and column names are the FITS keywords for
-    that metadata attribute.  For example, values from `meta.observation.time`
-    would be stored in the `TIME-OBS` column.
-
-    Rules for what function to use to determine the blended output attribute
-    value and what metadata attributes should be used as columns in the
-    blended FITS table are defined in the datamodel schema.
-
-
-    NOTE
-    ====
-    Custom rules for a metadata value should be computed by the calling routine
-    and used to update the metadata in the output model AFTER
-    calling this function.
-
-    Parameters
-    ----------
-    product : str
-        Name of combined product with metadata that needs updating. This can
-        be specified as a single filename.
-        When no value for `inputs` has been provided, this file
-        will also evaluate `meta.asn` to determine the names of the input
-        datamodels whose metadata need to be blended to create
-        the new combined metadata.
-
-    inputs : list, optional
-        This can be either a list of filenames or a list of DataModels objects.
-        If provided, the filenames provided in this list will be used to get
-        the metadata which will be blended into the final output metadata.
-
-    output : str, optional
-        If provided, update `meta.filename` in the blended `product`
-        to define what file this model will get written out to.
-
-    ignore : list of str, None, optional
-        A list of string the meta attribute names which, if provided,
-        will show which attributes should not be blended.
-
-    verbose : bool, optional [Default: False]
-        Print out additional messages during processing when specified.
-
-    Example
-    -------
-    This example shows how to blend the metadata from a set of DataModels
-    already read in memory for the product created by the `resample` step.
-    This example relies on the Association file used as the input to the
-    `resample` step to specify all the inputs for blending using the
-    following syntax::
-
-    >>> from stdatamodels.jwst import datamodels
-    >>> asnfile = "jw99999-a3001_20170327t121212_coron3_001_asn.json"
-    >>> data = datamodels.open(asnfile)
-    >>> input_models = [data[3], data[4]]  # we know the last datasets are SCIENCE
-    >>> blendmodels(data.meta.asn_table.products[0].name, inputs=input_models)
-    """
-    if inputs in EMPTY_LIST:
-        input_filenames = extract_filenames_from_product(product)
-        inputs = [datamodels.open(i) for i in inputs]  # return datamodels
-    else:
-        if isinstance(inputs, datamodels.JwstDataModel):
-            input_filenames = [i.meta.filename for i in inputs]
-        else:
-            input_filenames = inputs  # assume list of filenames as input
-    if verbose:
-        print('Creating blended metadata from: ')
-        for i in input_filenames:
-            print('\t{}'.format(i))
-
-    newmeta, newtab = get_blended_metadata(inputs, verbose=verbose)
-
-    # open product to update metadata from input models used to create product
-    if isinstance(product, str):
-        output_model = datamodels.open(product)
-    else:
-        output_model = product
-
-    '''
-    NOTE 17-Jan-2017:
-     Effort needs to be made to insure that new blended values conform
-     to the definitions of the attribute as provided by the schema.
-     This will address #1650 for the jwst package.
-       https://github.com/STScI-JWST/jwst/issues/1650
-    '''
-
-    # Start by identifying elements of the model which need to be ignored
-    ignore_list = _build_schema_ignore_list(newmeta._schema)
-    ignore_list += ['meta.wcs']  # Necessary since meta.wcs is not in schema
-    if ignore:
-        ignore_list.extend(ignore)
-
-    # Now assign values from new_hdrs to output_model.meta
-    flat_new_metadata = newmeta.to_flat_dict()
-
-    for attr in flat_new_metadata:
-        attr_use = not [attr.startswith(i) for i in ignore_list].count(True)
-        if attr.startswith('meta') and attr_use:
-            try:
-                output_model[attr] = newmeta[attr]
-            except KeyError:
-                # Ignore keys that are in the asdf tree but not in the schema
-                pass
-
-    # Apply any user-specified filename for output product
-    if output:
-        output_model.meta.filename = output
-    else:
-        # Otherwise, determine output filename from metadata
-        output = output_model.meta.filename
-
-    # Now, append HDRTAB as new element in datamodel
-    newtab_schema = build_tab_schema(newtab)
-    if hasattr(output_model, 'hdrtab'):
-        del output_model.hdrtab
-    output_model.add_schema_entry('hdrtab', newtab_schema)
-    output_model.hdrtab = fits_support.from_fits_hdu(newtab, newtab_schema)
-
-    # Clean up for the next run
-    del newmeta, newtab
-
-
-def get_blended_metadata(input_models, verbose=False):
-    """
-    Return a blended metadata instance and table based on the input datamodels.
-    This will serve as the primary interface for blending datamodels.
-
-    Parameters
-    ----------
-    input_models : list
-        Either a single list of filenames from which to extract the metadata to
-        be blended, or a list of `datamodels.JwstDataModel` objects to be blended.
-        The input models are assumed to have the blending rules defined as
-        an integral part of the schema definition for the model.
-
-    Returns
-    -------
-    metadata : list
-        A list of blended metadata instances, one for each i
-
-    new_table : object
-        Single fits.TableHDU object that contains the combined results from
-        all input headers(extension). Each row will correspond to an image,
-        and each column corresponds to a single keyword listed in the rules.
-
-    """
-    if not isinstance(input_models, list) and \
-       not isinstance(input_models, ModelContainer):
-        input_models = [input_models]
-
-    # Turn input filenames into a set of metadata objects
-    close_models = False
-    if isinstance(input_models[0], str):
-        # convert `input_models` to a list of datamodels
-        input_models = [datamodels.open(i) for i in input_models]
-        close_models = True
-
-    num_files = len(input_models)
-
-    # Determine what blending rules need to be merged to create the final
-    # blended headers. There will be a separate set of rules for each
-    # instrument, and all rules get merged into a composite set of rules that
-    # get applied to all input headers regardless of instrument.
-    #
-    # Instrument identification will be extracted from the INSTRUME keyword
-    # from the PRIMARY header of each input
-    #
-    icache = {}
-    for i in range(num_files):
-        model = input_models[i]
-        inst = model.meta.instrument.name.lower()
-
-        if inst not in icache:
-            # initialize the appropriate class for this data's instrument
-            inst_class = KeywordRules(model)
-            # Interpret rules for this class based on image that
-            # initialized this instrument's rules
-            inst_class.interpret_rules(model)
-
-            # Now add this interpreted class to the cache
-            icache[inst] = inst_class
-
-    # Create final merged set of rules
-    final_rules = None
-    for inst in icache:
-        if final_rules is None:
-            final_rules = icache[inst]
-        else:
-            final_rules.merge(icache[inst])
-
-    # apply rules to datamodels metadata
-    new_meta, newtab = final_rules.apply(input_models)
-
-    if len(newtab) > 0:
-        # Now merge the results for all the tables into single table extension
-        new_table = fits.BinTableHDU.from_columns(newtab)
-        new_table.header['EXTNAME'] = 'HDRTAB'
-    else:
-        new_table = None
-
-    if close_models:
-        for model in input_models:
-            model.close()
-
-    return new_meta, new_table
-
-
-def cat_headers(hdr1, hdr2):
+def blendmodels(product, inputs, ignore=None):
     """
-    Create new `astropy.io.fits.Header` object from
-    concatenating 2 input Headers
-    """
-    nhdr = hdr1.copy()
-    for c in hdr2.cards:
-        nhdr.append(c)
-
-    return fits.Header(nhdr)
-
-
-def extract_filenames_from_product(product):
-    """
-    Returns the list of filenames with extensions of input observations that
-    were used to generate the product.
-    """
-    asn_table = product.meta.asn.table_name
-    asn = associations.load_asn(asn_table)
-    prod = asn['products'][0]
-    fnames = [m['expname'] for m in prod['members']
-              if m['exptype'].upper() == 'SCIENCE']
-    return fnames
-
-
-def build_tab_schema(new_table):
-    """Return new schema definition that describes the input table."""
-    hdrtab = OrderedDict()
-    hdrtab['title'] = 'Combined header table'
-    hdrtab['fits_hdu'] = 'HDRTAB'
-    datatype = []
-    for col in new_table.columns:
-        cname = col.name
-        ctype = convert_dtype(str(col.dtype))
-        c = OrderedDict()
-        c['name'] = cname
-        c['datatype'] = ctype
-        datatype.append(c)
-    hdrtab['datatype'] = datatype
-
-    return hdrtab
-
-
-def convert_dtype(value):
-    """Convert numarray column dtype into YAML-compatible format description"""
-    if 'S' in value:
-        # working with a string description
-        str_len = int(value[value.find('S') + 1:])
-        new_dtype = ['ascii', str_len]
-    else:
-        new_dtype = str(value)
-
-    return new_dtype
-
-
-def _build_schema_ignore_list(schema):
-    """ Create a list of metadata that should be ignored when blending.
+    Blend datamodel metadata.
 
     Parameters
     ----------
-    schema : JSON schema fragment
-        The schema in which to search.
-
-    Returns
-    -------
-    results : list
-        List with schema attributes that needs to be ignored
-
-    """
-    def build_rules_list(subschema, path, combiner, ctx, recurse):
-        # Only interpret elements of the meta component of the model
-        if len(path) > 1 and path[0] == 'meta' and 'items' not in path:
-            attr = '.'.join(path)
-            if subschema.get('properties'):
-                return  # Ignore ObjectNodes
-            kwtype = subschema.get('type')
-            if kwtype == 'array':
-                results.append(attr)
-        else:
-            return
-
-    results = []
-    dm_schema.walk_schema(schema, build_rules_list, results)
-    return results
+    product: JwstDataModel
+        A datamodel that will have its metadata set
+        to the blended metadata and have the metadata
+        table assigned to the "hdrtab" attribute.
+
+    inputs: list of JwstDataModel
+        Input datamodels with metadata to blend.
+
+    ignore: list of str, optional
+        A list of metadata attributes to ignore during blending.
+        These attributes will not be set on the output/combined.
+        These attributes must be strings containing the dotted
+        path of each attribute (for example "meta.filename").
+        (Note that "meta.wcs" will always be ignored).
+    """
+    blender = ModelBlender(blend_ignore_attrs=ignore)
+    for model in inputs:
+        blender.accumulate(model)
+    blender.finalize_model(product)
diff --git a/jwst/model_blender/blendrules.py b/jwst/model_blender/blendrules.py
deleted file mode 100644
index 32bfd570a..000000000
--- a/jwst/model_blender/blendrules.py
+++ /dev/null
@@ -1,541 +0,0 @@
-""" blendmeta - Merge metadata from multiple models to create
-                   a new metadata instance and table
-
-"""
-from collections import OrderedDict
-
-import numpy as np
-from datetime import time
-from astropy.io import fits
-from astropy.time import Time
-from stdatamodels import schema as dm_schema
-from stdatamodels import properties
-
-from jwst import __version__
-from . import blender
-
-# Version of rules file format supported by this version of the code
-# All changes should be backwards compatible to older rules versions
-# so any rules file with Version >= __rules_version__ should work
-# with this code
-__rules_version__ = 2.1
-
-__all__ = ['find_keywords_in_section', 'first', 'float_one', 'int_one',
-           'interpret_attr_line', 'interpret_entry', 'last', 'maxdate',
-           'maxdatetime', 'maxtime', 'mindate', 'mindatetime', 'mintime',
-           'multi', 'multi1', 'zero', 'KeywordRules', 'KwRule']
-
-
-# Custom blending functions
-def multi(vals):
-    """
-    This will either return the common value from a list of identical values
-    or 'MULTIPLE'
-    """
-    uniq_vals = list(set(vals))
-    num_vals = len(uniq_vals)
-    if num_vals == 0:
-        return None
-    if num_vals == 1:
-        return uniq_vals[0]
-    if num_vals > 1:
-        return "MULTIPLE"
-
-
-def multi1(vals):
-    """
-    This will either return the common value from a list of identical values
-    or the single character '?'
-    """
-    uniq_vals = list(set(vals))
-    num_vals = len(uniq_vals)
-    if num_vals == 0:
-        return None
-    if num_vals == 1:
-        return uniq_vals[0]
-    if num_vals > 1:
-        return "?"
-
-
-def float_one(vals):
-    """ Return a constant floating point value of 1.0
-    """
-    return 1.0
-
-
-def int_one(vals):
-    """ Return an integer value of 1
-    """
-    return int(1)
-
-
-def zero(vals):
-    """ Return a value of 0
-    """
-    return 0
-
-
-def first(items):
-    """ Return first item from list of values"""
-    if len(items):
-        return items[0]
-    return None
-
-
-def last(items):
-    """ Return last item from list of values"""
-    if len(items):
-        return items[-1]
-    return None
-
-
-def mindate(items):
-    """Return the minimum date from a list of date strings in yyyy-mm-dd format."""
-    time_list = Time(items, format="iso", in_subfmt="date", out_subfmt="date")
-    return str(time_list.min())
-
-
-def maxdate(items):
-    """Return the maximum date from a list of date strings in yyyy-mm-dd format."""
-    time_list = Time(items, format="iso", in_subfmt="date", out_subfmt="date")
-    return str(time_list.max())
-
-
-def mindatetime(items):
-    """Return the minimum datetime from a list of datetime strings in ISO-8601 format."""
-    time_list = Time(items, format="isot")
-    return str(time_list.min())
-
-
-def maxdatetime(items):
-    """Return the maximum datetime from a list of datetime strings in ISO-8601 format."""
-    time_list = Time(items, format="isot")
-    return str(time_list.max())
-
-
-def mintime(items):
-    times = [_isotime(time_str) for time_str in items]
-    return min(times).isoformat()
-
-
-def maxtime(items):
-    times = [_isotime(time_str) for time_str in items]
-    return max(times).isoformat()
-
-
-def _isotime(time_str):
-    hms = [float(i) for i in time_str.split(':')]
-    sec_ms = hms[2] - int(hms[2])
-    isotime = time(int(hms[0]), int(hms[1]), int(hms[2]), int(sec_ms * 1000000))
-    return isotime
-
-
-# translation dictionary for function entries from rules files
-blender_funcs = {'first': first,
-                 'last': last,
-                 'float_one': float_one,
-                 'int_one': int_one,
-                 'zero': zero,
-                 'multi': multi,
-                 'multi?': multi1,
-                 'mean': np.mean,
-                 'sum': np.sum,
-                 'max': np.max,
-                 'min': np.min,
-                 'stddev': np.std,
-                 'mintime': mintime,
-                 'maxtime': maxtime,
-                 'mindate': mindate,
-                 'maxdate': maxdate,
-                 'mindatetime': mindatetime,
-                 'maxdatetime': maxdatetime}
-
-
-# Classes for managing keyword rules
-class KeywordRules():
-
-    def __init__(self, model):
-        """ Read in the rules used to interpret the keywords from the specified
-            instrument image header.
-        """
-        self.instrument = model.meta.instrument.name.lower()
-        self.new_header = None
-        self.rules_version = __rules_version__
-
-        self.rules_file = model._schema
-        self.rule_specs = _build_schema_rules_dict(self.rules_file)
-
-        self.rule_objects = []
-        self.rules = []
-        self.section_names = []
-
-    def interpret_rules(self, hdrs):
-        """ Convert specifications for rules from rules file
-        into specific rules for this header(instrument/detector).
-
-        Notes
-        -----
-        This allows for expansion rules to be applied to rules
-        from the rules files (such as any wildcards or section titles).
-
-        Output will be 'self.rules' that contains a list of tuples:
-        - a tuple of 2 values for each column in the table
-        - a tuple of 4 values for each attribute identified in metadata
-        Partial sample from HST to show format:
-        [('CTYPE1O', 'CTYPE1O'),
-        ('CTYPE2O', 'CTYPE2O'),
-        ('CUNIT1O', 'CUNIT1O'),
-        ('CUNIT2O', 'CUNIT2O'),
-        ('APERTURE', 'APERTURE', <function fitsblender.blendheaders.multi>, 'ignore'),
-        ('DETECTOR', 'DETECTOR', <function fitsblender.blender.first>, 'ignore'),
-        ('EXPEND', 'EXPEND', <function numpy.core.fromnumeric.amax>, 'ignore'),
-        ('EXPSTART', 'EXPSTART', <function numpy.core.fromnumeric.amin>, 'ignore'),
-        ('EXPTIME', 'TEXPTIME', <function numpy.core.fromnumeric.sum>, 'ignore'),
-        ('EXPTIME', 'EXPTIME', <function numpy.core.fromnumeric.sum>, 'ignore')]
-
-        This rules format will allow the algorithm, logic and code from
-        the original fitsblender to be used with as little change as
-        possible.  It will need to be derived (as with HST) from the
-        input models metadata for expansion of attribute sections or
-        wildcards in attributes specified in the rules.
-
-        """
-        if isinstance(hdrs, tuple):
-            hdrs = list(hdrs)
-        if not isinstance(hdrs, list):
-            hdrs = [hdrs]
-
-        # apply rules to headers
-        for attr in self.rule_specs:
-            speclist = self.rule_specs[attr]
-
-            for rule in speclist:
-                # Create KwRule input equivalent to HST rules input
-                kwr = KwRule(rule)
-                duplicate_rule = False
-
-                for robj in self.rule_objects:
-                    if kwr.rule_spec == robj.rule_spec:
-                        duplicate_rule = True
-                        break
-                if not duplicate_rule:
-                    for hdr in hdrs:
-                        kwr.interpret(hdr)
-
-                    self.rules.extend(kwr.rules)
-
-    def merge(self, kwrules):
-        """
-        Merge a new set of interpreted rules into the current set
-        The new rules, kwrules, can either be a new class or a whole new
-        set of rules (like those obtained from using self.interpret_rules with
-        a new header).
-        """
-        if isinstance(kwrules, KeywordRules):
-            kwrules = kwrules.rules
-
-        # Determine what rules are specified in kwrules that
-        #    are NOT in self.rules
-        k = []
-        # Delete these extraneous rules from input kwrules
-        for r in kwrules:
-            if r not in self.rules:
-                k.append(r)
-
-        # extend self.rules with additional rules
-        self.rules.extend(k)
-
-    def apply(self, models, tabhdu=False):
-        """ For a full list of metadata objects, apply the specified rules to
-            generate a dictionary of new values and a table using
-            blender.
-
-            This method returns the new metadata object and summary table
-            as `datamodels.model.ndmodel` and fits.binTableHDU objects.
-        """
-        # Apply rules to headers
-        fbdict, fbtab = blender.metablender(models, self.rules)
-
-        # Determine which keywords are included in the table but not
-        # the new dict(header). These will be removed from the output
-        # header altogether
-        tabcols = fbtab.dtype.names
-
-        # Start with a copy of the template as the new header
-        # This will define what keywords need to be updated, as the rules
-        # and input headers often include headers for multiple extensions in
-        # order to build the complete table for all the keywords in the file
-        # in one run
-        new_model = models[0].copy()
-
-        # Apply updated/blended values into new header, but only those
-        # keywords which are already present in the 'template' new header
-        # this allows the rules to be used on all extensions at once yet
-        # update each extension separately without making copies of kws from
-        # one extension to another.
-        # new_model.update(fbdict)
-        for attr in new_model.to_flat_dict():
-            if 'meta' in attr and attr in fbdict:
-                new_model[attr] = fbdict[attr]
-
-        # Create summary table
-        if len(tabcols) > 0:
-            if tabhdu:
-                new_table = fits.BinTableHDU.from_columns(fbtab)
-                new_table.header['EXTNAME'] = 'HDRTAB'
-            else:
-                new_table = fbtab
-        else:
-            new_table = None
-        return new_model, new_table
-
-    def add_rules_kws(self, hdr):
-        """ Update metadata with
-        WARNING
-        -------
-        Needs to be modified to work with metadata.
-
-        Update PRIMARY header with HISTORY cards that report the exact
-        rules used to create this header. Only non-comment lines from the
-        rules file will be reported.
-        """
-        hdr['meta.rules_version'] = (__rules_version__,
-                                     'Version ID for header kw rules file')
-        hdr['meta.blender_version'] = (__version__,
-                                       'Version of blendheader software used')
-
-    def index_of(self, kw):
-        """ Reports the index of the specified kw."""
-        indx = []
-        for r, i in zip(self.rules, list(range(len(self.rules)))):
-            if r[0] == kw:
-                indx.append(i)
-        return indx
-
-
-class KwRule():
-    """
-    This class encapsulates the logic needed for interpreting a single keyword
-    rule from a text file.
-
-    Notes
-    -----
-    The ``.rules`` attribute contains the interpreted set of rules that corresponds
-    to this line.
-
-    Example::
-
-      Interpreting rule from
-      {'meta.attribute': { 'rule': 'first', 'output': 'meta.attribute'}}
-      --or--
-      {'meta.attribute': 'meta.attribute'}  # Table column specification
-
-      into rule [('meta.attribute', 'meta.attribute', <function first at 0x7fe505db7668>, 'ignore')]
-      and sname None
-
-    """
-
-    def __init__(self, line):
-        """Initialize new keyword rule.
-
-        Parameters
-        ==========
-        line : dict
-            Line should be dict with attribute name as the key, and
-            a dict as the value specifying 'rule' and (optionally)'output'.
-        """
-        self.rule_spec = line  # dict read in from rules file
-        self.rules = []
-        self.delete_kws = []
-        self.section_name = []
-
-    def interpret(self, hdr):
-        """Use metadata to interpret rule."""
-        if self.rules:
-            # If self.rules has already been defined for this rule, do not try
-            # to interpret it any further with additional headers
-            return
-        irules, sname = interpret_entry(self.rule_spec, hdr)
-
-        # keep track of any section name identified for this rule
-        if sname:
-            self.section_name.append(sname)
-
-        # Now, interpret rule based on presence of kw in hdr
-        if irules:
-            self.rules = irules
-
-
-# Utility functions.
-def _build_schema_rules_dict(schema):
-    """ Create a dict that extracts blend rules from an input schema.
-
-    Parameters
-    ----------
-    schema : JSON schema fragment
-        The schema in which to search.
-
-    Returns
-    -------
-    results : OrderedDict
-        Dictionary with schema attributes as keys and blend rules
-        as values
-
-    """
-    def build_rules_dict(subschema, path, combiner, ctx, recurse):
-        # Only interpret elements of the meta component of the model
-        if len(path) > 1 and path[0] == 'meta' and 'items' not in path:
-            for combiner in ['anyOf', 'oneOf']:
-                if combiner in path:
-                    path = path[:path.index(combiner)]
-                    break
-            attr = '.'.join(path)
-            if subschema.get('properties'):
-                return  # Ignore ObjectNodes
-
-            # Get blending info
-            kwrule = subschema.get('blend_rule')
-            kwtab = subschema.get('blend_table')
-            kwname = subschema.get('fits_keyword', attr)
-
-            # If rules had already been set, only modify
-            # the rules if there are explicit settings.
-            rule_spec = None
-            result = results.get(attr, [])
-            if kwrule:
-                rule_spec = {attr: {'rule': kwrule}}
-            elif not results.get(attr):
-                rule_spec = {attr: {'rule': 'first'}}
-            if rule_spec:
-                result.append(rule_spec)
-
-            # Add keyword to table if specified.
-            if kwtab:
-                result.append({attr: kwname})
-
-            # Add the results back.
-            if len(result):
-                results[attr] = result
-
-        else:
-            return
-
-    results = OrderedDict()
-    dm_schema.walk_schema(schema, build_rules_dict, results)
-    return results
-
-
-def interpret_entry(line, hdr):
-    """ Generate the rule(s) specified by the entry from the rules file.
-
-    Notes
-    -----
-    The entry should always be a dict with format:
-    {attribute_name : {'rule':'some_rule', 'output':''}}
-    -- or (for table column specification)--
-    {attribute_name: attribute_name}
-    where
-    'output' is assumed to be the same as attribute_name if not present
-
-    """
-    # Interpret raw input line
-    attr = list(line.keys())[0]
-    line_spec = line[attr]
-    attr_column = True  # Determine whether this rule defines a table column
-    if isinstance(line_spec, dict):
-        attr_column = False  # If not, turn this off
-
-    # Initialize output values
-    rules = []
-    section_name = None
-
-    # Parse the line
-    if isinstance(hdr[attr], properties.ObjectNode):
-        section_name = attr
-        # Datamodel sections are just parent Nodes for each attribute
-        keys = hdr[section_name].instance.keys()
-        kws = ['{}.{}'.format(section_name, k) for k in keys]
-        #
-        if kws is not None:
-            for kw1 in kws:
-                if attr_column:
-                    # Interpret each attribute as a new column in Table
-                    rules.append((kw1, kw1))
-                else:
-                    # Interpret each attribute in Node with same rule
-                    attr_rules = interpret_attr_line(kw1, line_spec)
-                    rules.extend(attr_rules)
-    else:
-        attr_rules = interpret_attr_line(attr, line_spec)
-        rules.extend(attr_rules)
-
-    return rules, section_name
-
-
-def interpret_attr_line(attr, line_spec):
-    """ Generate rule for single attribute from input line from rules file."""
-    rules = []
-
-    kws = [attr]
-    if isinstance(line_spec, dict):
-        if 'output' in line_spec:
-            kws2 = [line_spec['output']]
-        else:
-            kws2 = kws
-    else:
-        kws2 = [line_spec]
-
-    lrule = None
-    if 'rule' in line_spec:
-        lrule = line_spec['rule']
-
-    # Interpret short-hand rules using dict
-    if lrule is not None and len(lrule) > 0:
-        if lrule in blender_funcs:
-            lrule = blender_funcs[lrule]
-        else:
-            lrule = None
-        # build separate rule for each kw
-        for kw1, kw2 in zip(kws, kws2):
-            new_rule = (kw1, kw2, lrule, "ignore")
-            if new_rule not in rules:
-                rules.append(new_rule)
-    else:
-        for kw1, kw2 in zip(kws, kws2):
-            new_rule = (kw1, kw2)
-            if new_rule not in rules:
-                rules.append(new_rule)
-
-    return rules
-
-
-def find_keywords_in_section(hdr, title):
-    """Return a list of keyword names.
-
-    The list will be derived from the section
-        with the specified section title identified in the hdr.
-    """
-    # Identify card indices of start and end of specified section
-    sect_start = None
-    sect_end = None
-    for i, kw in enumerate(hdr.cards):
-        if sect_start is None:
-            if title in str(hdr[i]):
-                sect_start = i
-        else:
-            if '/' in str(hdr[i]) and hdr[i] not in ['N/A', ' ', '']:
-                sect_end = i
-                break
-    if sect_end is None:
-        sect_end = len(hdr)
-    if sect_start is None:
-        return None
-
-    # Now, extract the keyword names from this section
-    # section_keys = hdr.ascard[sect_start+1:sect_end-1].keys()
-    section_keys = list(hdr[sect_start + 1:sect_end - 1].keys())
-    # remove any blank keywords
-    while section_keys.count('') > 0:
-        section_keys.remove('')
-
-    return section_keys
diff --git a/jwst/model_blender/rules.py b/jwst/model_blender/rules.py
new file mode 100644
index 000000000..799d6ffba
--- /dev/null
+++ b/jwst/model_blender/rules.py
@@ -0,0 +1,112 @@
+import numpy as np
+
+
+__all__ = ['RULE_FUNCTIONS', 'AttributeBlender', 'make_blender']
+
+
+def _multi(vals):
+    """
+    This will either return the common value from a list of identical values
+    or 'MULTIPLE'
+    """
+    uniq_vals = list(set(vals))
+    num_vals = len(uniq_vals)
+    if num_vals == 0:
+        return None
+    if num_vals == 1:
+        return uniq_vals[0]
+    if num_vals > 1:
+        return "MULTIPLE"
+
+
+RULE_FUNCTIONS = {
+    'multi': _multi,
+    'mean': np.mean,
+    'sum': np.sum,
+    'max': np.max,
+    'min': np.min,
+
+     # retained date/time names for backwards compatibility
+     # as these all assume ISO8601 format the lexical and
+     # chronological sorting match
+    'mintime': min,
+    'maxtime': max,
+    'mindate': min,
+    'maxdate': max,
+    'mindatetime': min,
+    'maxdatetime': max
+}
+"""
+Mapping of rule names to functions.
+
+Used for `make_blender`.
+
+The following rules are considered deprecated
+and should not be used for new schemas.
+
+  - mintime
+  - maxtime
+  - mindate
+  - maxdate
+  - mindatetime
+  - maxdatetime
+"""
+
+
+class AttributeBlender:
+    """
+    Single attribute metadata blender
+    """
+    def __init__(self, blend_function):
+        """
+        Create a new metadata attribute blender.
+
+        Parameters
+        ----------
+        blend_function: callable
+            Function to blend accumulated metadata values
+        """
+        self.blend_function = blend_function
+        self.values = []
+
+    def accumulate(self, value):
+        """
+        Add a metadata value for blending.
+
+        Parameters
+        ----------
+        value:
+            Value for this metadata attribute to use
+            when blending.
+        """
+        self.values.append(value)
+
+    def finalize(self):
+        """
+        Blend the accumulated metadata values.
+
+        Returns
+        -------
+        value:
+            The blended result.
+        """
+        if not self.values:
+            return None
+        return self.blend_function(self.values)
+
+
+def make_blender(rule):
+    """
+    Make a `AttributeBlender` instance using the provided rule
+
+    Parameters
+    ----------
+    rule: string
+        Name of the blending rule. Must be in `RULE_FUNCTIONS`.
+
+    Returns
+    -------
+    attr_blender: `AttrBlender`
+        Blender instance using the provided rule.
+    """
+    return AttributeBlender(RULE_FUNCTIONS[rule])
diff --git a/jwst/model_blender/textutil.py b/jwst/model_blender/textutil.py
deleted file mode 100644
index b5b14a541..000000000
--- a/jwst/model_blender/textutil.py
+++ /dev/null
@@ -1,72 +0,0 @@
-"""Text output-related utilities."""
-
-
-import textwrap
-
-
-def wrap(text, width, *args, **kwargs):
-    """
-    Like :func:`textwrap.wrap` but preserves existing newlines which
-    :func:`textwrap.wrap` does not otherwise handle well.
-
-    See Also
-    --------
-    :func:`textwrap.wrap`
-    """
-
-    return sum([textwrap.wrap(line, width, *args, **kwargs)
-                if line else [''] for line in text.splitlines()], [])
-
-
-def textbox(text, width=78, boxchar='#', indent=0):
-    """
-    Outputs line-wrapped text wrapped in a box drawn with a repeated (usually
-    ASCII) character.
-
-    For example:
-
-        >>> print(textbox('Text to wrap', width=16))
-        ################
-        #              #
-        # Text to wrap #
-        #              #
-        ################
-
-    Parameters
-    ----------
-    text : string
-        The text to wrap
-
-    width : int
-        The width of the entire box, including the perimeter and
-        the indentation space.  Because the
-        wrapped text is padded with an additional column of whitespace on each
-        side, the minimum width is 5--any width less than that is
-        is automatically increased to 5 (default: 78)
-
-    boxchar : string
-        (No pun intended.) The character to draw the box with.  May also
-        be a string of multiple characters (default: '#')
-
-    indent : int
-        Amount of space by which the box should be indented. (default: 0)
-
-    """
-
-    min_width = len(boxchar) * 2 + 3
-    width = max(width - indent, min_width)
-    indentspace = indent * ' '
-
-    wrap_width = width - min_width + 1
-
-    q, r = divmod(width, len(boxchar))
-    # The top/bottom border
-    top_border = indentspace + boxchar * q + boxchar[:r]
-    top_padding = indentspace + boxchar + ' ' * (width - len(boxchar) * 2) + boxchar
-
-    lines = ['%s%s %s %s' % (indentspace, boxchar, line.ljust(wrap_width),
-                             boxchar)
-             for line in wrap(text, wrap_width)]
-    top = [top_border, top_padding]
-    bottom = [top_padding, top_border]
-    return '\n'.join(top + lines + bottom)
diff --git a/jwst/pipeline/calwebb_coron3.py b/jwst/pipeline/calwebb_coron3.py
index 04673418a..b35cceecf 100644
--- a/jwst/pipeline/calwebb_coron3.py
+++ b/jwst/pipeline/calwebb_coron3.py
@@ -7,7 +7,7 @@ from stdatamodels.jwst import datamodels
 
 from jwst.datamodels import ModelContainer, ModelLibrary
 
-from ..model_blender import blendmeta
+from ..model_blender import ModelBlender
 
 # step imports
 from ..coron import stack_refs_step
@@ -161,8 +161,10 @@ class Coron3Pipeline(Pipeline):
         # Call the sequence of steps: outlier_detection, align_refs, and klip
         # once for each input target exposure
         resample_input = ModelContainer()
+        model_blender = ModelBlender()
         for target_file in targ_files:
             with datamodels.open(target_file) as target:
+                model_blender.accumulate(target)
 
                 # Remove outliers from the target
                 if not skip_outlier_detection:
@@ -213,7 +215,7 @@ class Coron3Pipeline(Pipeline):
             completed = 'SKIPPED'
         if completed == 'COMPLETE':
             self.log.debug(f'Blending metadata for {result}')
-            blendmeta.blendmodels(result, inputs=targ_files)
+            model_blender.finalize_model(result)
 
         try:
             result.meta.asn.pool_name = input_models.meta.asn_table.asn_pool
diff --git a/jwst/resample/resample.py b/jwst/resample/resample.py
index 97b2e3d76..6350c3496 100644
--- a/jwst/resample/resample.py
+++ b/jwst/resample/resample.py
@@ -15,8 +15,8 @@ from jwst.datamodels import ModelLibrary
 from jwst.associations.asn_from_list import asn_from_list
 
 from . import gwcs_drizzle
+from jwst.model_blender.blender import ModelBlender
 from jwst.resample import resample_utils
-from ..model_blender import blendmeta
 
 log = logging.getLogger(__name__)
 log.setLevel(logging.DEBUG)
@@ -201,24 +201,6 @@ class ResampleData:
         else:
             return self.resample_many_to_one(input_models)
 
-    def blend_output_metadata(self, output_model, input_models):
-        """Create new output metadata based on blending all input metadata."""
-        # Run fitsblender on output product
-        output_file = output_model.meta.filename
-
-        ignore_list = [
-            'meta.photometry.pixelarea_steradians',
-            'meta.photometry.pixelarea_arcsecsq',
-        ]
-
-        log.info('Blending metadata for {}'.format(output_file))
-        blendmeta.blendmodels(
-            output_model,
-            inputs=input_models,
-            output=output_file,
-            ignore=ignore_list
-        )
-
     def _get_intensity_scale(self, img):
         """
         Compute an intensity scale from the input and output pixel area.
@@ -387,28 +369,13 @@ class ResampleData:
         copy_asn_info_from_library(input_models, output_model)
 
         if self.blendheaders:
-            # right now this needs a list of input models, all in memory
-            # for now, just load the models as a list with empty data arrays
-            # but the blend_meta step itself should eventually be refactored 
-            # to expect a list of metadata objects
-            # instead of a list of datamodels
-            input_list = []
-            with input_models:
-                for i, model in enumerate(input_models):
-                    empty_model = type(model)()
-                    empty_model.meta = model.meta
-                    copy_asn_info_from_library(input_models, empty_model)
-                    empty_model.data = np.empty((1, 1))
-                    empty_model.dq = np.empty((1, 1))
-                    empty_model.err = np.empty((1, 1))
-                    empty_model.wht = np.empty((1, 1))
-                    empty_model.var_rnoise = np.empty((1, 1))
-                    empty_model.var_poisson = np.empty((1, 1))
-                    empty_model.var_flat = np.empty((1, 1))
-                    input_list.append(empty_model)
-                    input_models.shelve(model, i, modify=False)
-            self.blend_output_metadata(output_model, input_list)
-            del input_list
+            blender = ModelBlender(
+                blend_ignore_attrs=[
+                    'meta.photometry.pixelarea_steradians',
+                    'meta.photometry.pixelarea_arcsecsq',
+                    'meta.filename',
+                ]
+            )
 
         # Initialize the output with the wcs
         driz = gwcs_drizzle.GWCSDrizzle(output_model, pixfrac=self.pixfrac,
@@ -417,6 +384,8 @@ class ResampleData:
         log.info("Resampling science data")
         with input_models:
             for img in input_models:
+                if self.blendheaders:
+                    blender.accumulate(img)
                 iscale = self._get_intensity_scale(img)
                 log.debug(f'Using intensity scale iscale={iscale}')
                 img.meta.iscale = iscale
@@ -449,6 +418,9 @@ class ResampleData:
                 del data, inwht
                 input_models.shelve(img)
 
+        if self.blendheaders:
+            blender.finalize_model(output_model)
+
         # Resample variance arrays in input_models to output_model
         self.resample_variance_arrays(output_model, input_models)
         var_components = [
@@ -979,4 +951,4 @@ def copy_asn_info_from_library(library, output_model):
     if (
         asn_table_name := library.asn.get("table_name", None)
     ) is not None:
-        output_model.meta.asn.table_name = asn_table_name
\ No newline at end of file
+        output_model.meta.asn.table_name = asn_table_name


Test patch:
diff --git a/jwst/model_blender/tests/test_blend.py b/jwst/model_blender/tests/test_blend.py
index dc0e578b7..37c9d705d 100644
--- a/jwst/model_blender/tests/test_blend.py
+++ b/jwst/model_blender/tests/test_blend.py
@@ -1,13 +1,13 @@
 """Test blend_models"""
 import pytest
 
-from astropy.table import Table
 import numpy as np
 
 from stdatamodels.schema import walk_schema
-from stdatamodels.jwst.datamodels import ImageModel
+from stdatamodels.jwst.datamodels import ImageModel, CubeModel
 
 from jwst.model_blender import blendmeta
+from jwst.model_blender.blender import ModelBlender
 
 
 # Setup various input meta data
@@ -19,8 +19,23 @@ FILENAMES = ['image1_cal.fits', 'image2_cal.fits', 'image3_cal.fits']
 DATETIMES = ['2017-11-30T13:52:20.367', '2017-11-11T15:14:29.176',
              '2017-11-11T15:15:06.118']
 DATES = ['2017-11-30', '2017-11-11', '2017-12-10']
-INSTRUMENT_NAMES = ['NIRCAM'] * 3
+INSTRUMENT_NAMES = ['NIRCAM'] * N_MODELS
+DETECTOR_NAMES = ['NRCA1', 'NRCA2', 'NRCA3']
 CORONAGRAPHS = ['4QPM', '4QPM_1065', '4QPM_1140']
+EXP_ONLYS = [True] * N_MODELS
+POLYNOMIAL_INFOS = [[{'degree': [1]}]] * N_MODELS
+
+# even though this won't be in the schema the blender
+# copies the metadata from the first model so the
+# meta from the first model will appear in the result
+UNKNOWN_IN_FIRSTS = ['foo', 'bar', 'bam']
+
+# but... if the metadata isn't in the first model it
+# won't appear in the output
+UNKNOWN_MISSING_FIRSTS = [None, 'fizz', 'buzz']
+
+# even if the metadata is defined in the schema
+KNOWN_MISSING_FIRSTS = [None, '1', '2']
 
 
 def _make_data():
@@ -46,20 +61,36 @@ def _make_data():
         'meta.filename': FILENAMES,
         'meta.instrument.coronagraph': CORONAGRAPHS,
         'meta.instrument.name': INSTRUMENT_NAMES,
+        'meta.instrument.detector': DETECTOR_NAMES,
         'meta.date': DATETIMES,
         'meta.observation.date': DATES,
         'meta.observation.date_beg': DATETIMES,
+        'meta.visit.exp_only': EXP_ONLYS,
+        'meta.background.polynomial_info': POLYNOMIAL_INFOS,
+
+        # test a few "quirks" of the blender including...
+        # unknown (not in schema) attribute in the first model
+        'meta.unknown_in_first': UNKNOWN_IN_FIRSTS,
+        # unknown attribute missing from the first model
+        'meta.unknown_missing_first': UNKNOWN_MISSING_FIRSTS,
+        # known (in schema) attribute missing from the first model
+        'meta.observation.visit_id': KNOWN_MISSING_FIRSTS,
     }
     output_values = {
         'meta.exposure.start_time': START_TIMES[0],
         'meta.exposure.exposure_time': np.sum(EXP_TIMES),
         'meta.exposure.end_time': END_TIMES[-1],
-        'meta.filename': FILENAMES[0],
+        # skip filename as it will be ignored below
         'meta.instrument.coronagraph': CORONAGRAPHS[0],
         'meta.instrument.name': INSTRUMENT_NAMES[0],
+        'meta.instrument.detector': 'MULTIPLE',
         'meta.date': DATETIMES[0],
         'meta.observation.date': DATES[1],
-        'meta.observation.date_beg': DATETIMES[1]
+        'meta.observation.date_beg': DATETIMES[1],
+        'meta.visit.exp_only': EXP_ONLYS[0],
+
+        # this is the only "quirk" that will appear in the output
+        'meta.unknown_in_first': UNKNOWN_IN_FIRSTS[0],
     }
 
     for i, model in enumerate(models):
@@ -88,8 +119,10 @@ def blend(make_data):
         or from the results of `_make_data` directly.
     """
     models, input_values, output_values = make_data
-    newmeta, newtab = blendmeta.get_blended_metadata(models)
-    return newmeta, newtab, models, input_values, output_values
+    output = ImageModel()
+    blendmeta.blendmodels(output, models, ignore=["meta.filename"])
+    newmeta = {k: v for k, v in output.to_flat_dict().items() if k.startswith('meta')}
+    return newmeta, output.hdrtab, models, input_values, output_values
 
 
 def test_blendmeta(blend):
@@ -103,7 +136,10 @@ def test_blendmeta(blend):
     newmeta, newtab, models, input_values, output_values = blend
 
     for attr in input_values:
-        assert newmeta[attr] == output_values[attr]
+        if attr in output_values:
+            assert newmeta[attr] == output_values[attr]
+        else:
+            assert attr not in newmeta
 
 
 def build_fits_dict(schema):
@@ -152,9 +188,20 @@ def test_blendtab(blend):
     meta_to_fits = dict(map(reversed, fits_to_meta.items()))
     fits_expected = set(
         meta_to_fits[meta]
-        for meta in input_values
+        for meta in output_values
+        if meta in meta_to_fits
     )
 
     # Ensure all the expected FITS keywords are in the table.
-    table = Table(newtab.data)
-    assert not fits_expected.difference(table.colnames)
+    colnames = set(newtab.dtype.fields)
+    assert not fits_expected.difference(colnames)
+    for col in colnames:
+        if col in input_values:
+            assert newtab[col] == input_values[col]
+
+
+def test_wrong_type_failure():
+    blender = ModelBlender()
+    blender.accumulate(ImageModel())
+    with pytest.raises(ValueError, match="model of type"):
+        blender.accumulate(CubeModel())