Map channels to Var Channels

src/opsinputs/VarObsWriter.cc

@@ -17,6 +17,7 @@
 #include "ioda/ObsVector.h"
 #include "oops/base/Variables.h"
 #include "oops/mpi/mpi.h"
+#include "oops/util/IntSetParser.h"
 #include "oops/util/Logger.h"
 #include "opsinputs/LocalEnvironment.h"
 #include "opsinputs/VarObsWriterParameters.h"
@@ -70,6 +71,13 @@ VarObsWriter::VarObsWriter(ioda::ObsSpace & obsdb, const Parameters_ & params,
   const int fallbackChannels = 0;
   // Avoid passing a null pointer to Fortran.
   const int *channelsData = channels.empty() ? &fallbackChannels : channels.data();
+
+  // Want to also set up an ordered list of channels numbered for VAR.
+  std::set<int> varchanset = oops::parseIntSet(parameters_.varChannels.value());
+  varchannels_.assign(varchanset.begin(), varchanset.end());
+  conf.set("varChannels", varchannels_);
+  conf.set("NumVarChannels", sizeof(varchannels_));
+
   if (!opsinputs_varobswriter_create_f90(key_, &conf,
                                          fortranMpiCommunicatorIsValid,
                                          fortranMpiCommunicator,

src/opsinputs/VarObsWriter.h

@@ -8,7 +8,9 @@
 
 #include <memory>
 #include <ostream>
+#include <set>
 #include <string>
+#include <vector>
 
 #include "ioda/ObsDataVector.h"
 #include "oops/base/Variables.h"
@@ -87,6 +89,7 @@ class VarObsWriter : public oops::interface::ObsFilterBase<ufo::ObsTraits>,
   oops::Variables extradiagvars_;
   std::shared_ptr<ioda::ObsDataVector<int>> flags_;
   std::shared_ptr<ioda::ObsDataVector<float>> obsErrors_;
+  std::vector<int> varchannels_;
 
   VarObsWriterParameters parameters_;
 };

src/opsinputs/VarObsWriterParameters.h

@@ -79,21 +79,25 @@ class VarObsWriterParameters : public oops::ObsFilterParametersBase {
   /// Update OPS flag to output the varbc predictors
   oops::Parameter<bool> outputVarBCPredictors{"output_varbc_predictors", false, this};
 
-  /// This contains the offset that needs to be added to the channel number in order to
-  /// index the output arrays correctly.
-  oops::Parameter<int> channel_offset{"channel_offset", 0, this};
-
   /// This is the size of the varobs array for output.  The default is zero and the size
   /// of the array will be used.
   /// For atovs, jedi has 20 brightness temperatures but var expects 40.
   /// Therefore for atovs brightness_tmperatuere => size_of_varobs_array = 40.
   oops::Parameter<int> size_of_varobs_array{"size_of_varobs_array", 0, this};
 
-  /// This matches the channel number with the array index. This is useful when
-  /// there are skipped channels. e.g. channels 5,6,7,9,10,11 are required from a possible 12
-  /// channels. This would fill a size 12 array with the array indexes matching the channel number
-  /// [NaN,NaN,NaN,Nan,5,6,7,Nan,9,10,11,NaN].
-  oops::Parameter<bool> use_actual_channels{"use_actual_channels", false, this};
+  /// List of channels which are expected for Var. This enables mapping of JOPA
+  /// channel numbers to VAR numbers. Default is empty and the JOPA channels will be used.
+  oops::Parameter<std::string> varChannels{"varChannels", "", this};
+
+  /// If set to false (default true) this matches the channel number with the array index.
+  /// This is useful when there are skipped channels. e.g. channels 5,6,7,9,10,11 are
+  /// required from a possible 12 channels. This would fill a size 12 array with the array
+  /// indexes matching the channel number [NaN,NaN,NaN,Nan,5,6,7,Nan,9,10,11,NaN].
+  /// Compressed example could be [2,8,50,100], no missing values between channel numbers.
+  oops::Parameter<bool> compressVarChannels{"compress_var_channels", true, this};
+
+  /// Increase the channel array size to the same size as the varobs array
+  oops::Parameter<bool> increaseChanArray{"increase_chan_array", false, this};
 
   /// If this list of ufo::variable is defined in the yaml a subset of the flags
   /// will be made with just these variables present.  This will allow Fortran calls such-as

src/opsinputs/opsinputs_fill_mod.F90

@@ -1030,12 +1030,19 @@ end subroutine opsinputs_fill_fillreal
 !>   of variables with suffixes corresponding to the indices specified in \p Channels.
 !> \param[in] JediGroup
 !>   Group of the JEDI variable used to populate \p Real2.
+!> \param[in] compressVarChannels
+!>   Whether to apply var channel compression (No NaN spaces between channels)
+!> \param[in] sizeOfVarobsArray
+!>   The size of the varobs array which the output data will be stored in.
+!> \param[in] varChannels
+!>   A list of the var channel numbers which the channels will be mapped to.
 !>
 !> \note This function returns early (without a warning) if the specified JEDI variable is not found.
 !> We rely on warnings printed by the OPS code whenever data needed to output a requested varfield
 !> are not found.
 subroutine opsinputs_fill_fillreal2d_norecords( &
-  Hdr, OpsVarName, NumObs, Real2, ObsSpace, Channels, JediVarName, JediVarGroup, OffsetChans, useActualChans)
+  Hdr, OpsVarName, NumObs, Real2, ObsSpace, Channels, JediVarName, &
+  JediVarGroup, compressVarChannels, sizeOfVarobsArray, varChannels)
 implicit none
 
 ! Subroutine arguments:
@@ -1047,8 +1054,10 @@ subroutine opsinputs_fill_fillreal2d_norecords( &
 integer(c_int), intent(in)                          :: Channels(:)
 character(len=*), intent(in)                        :: JediVarName
 character(len=*), intent(in)                        :: JediVarGroup
-type(opsinputs_channeloffset), optional, intent(in) :: OffsetChans
-logical, optional, intent(in)                       :: useActualChans
+logical, optional, intent(in)                       :: compressVarChannels
+integer(integer64), optional, intent(in)            :: sizeOfVarobsArray
+integer(c_int), optional, intent(in)                :: varChannels(:)
+
 
 ! Local declarations:
 real(kind=c_double)                             :: VarValue(NumObs)
@@ -1057,52 +1066,77 @@ subroutine opsinputs_fill_fillreal2d_norecords( &
 integer                                         :: iChannel
 integer                                         :: offset
 integer                                         :: numchans
-logical                                         :: localUseActualChans
+integer                                         :: offsetsize
+integer                                         :: arrayindex
+logical                                         :: compressChannels
 
 ! Body:
 
 MissingDouble = missing_value(0.0_c_double)
 
+compressChannels = .true.
+if (present(compressVarChannels)) then
+  compressChannels = compressVarChannels
+end if
+
 JediVarNamesWithChannels = opsinputs_fill_varnames_with_channels(JediVarName, Channels)
 
 !take into account offsetting of 2nd dimension if required
 !designed to be used to pack where multiple satellite instruments expected
 !e.g. HIRS in ATOVS stream
-offset = 0
-numchans = size(JediVarNamesWithChannels)
-if (present(OffsetChans)) then
-  offset = OffsetChans % channel_offset
-  if (OffsetChans % size_of_varobs_array > 0) &
-    numchans = OffsetChans % size_of_varobs_array
-end if
 
-!Setup for channels needing to match array index
-localUseActualChans = .false.
-if (present(useActualChans)) then
-    localUseActualChans = useActualChans
+numchans = size(JediVarNamesWithChannels)
+!sizeOfVarobsArray comes from intitial setting of size_of_varobs_array 
+! used to define the size of the channel array to fill. 
+if (present(sizeOfVarobsArray)) then
+  if (sizeOfVarobsArray > 0) then
+    numchans = sizeOfVarobsArray
+  end if
 end if
 
 if (obsspace_has(ObsSpace, JediVarGroup, JediVarNamesWithChannels(1))) then
   ! Allocate OPS data structures
   call Ops_Alloc(Hdr, OpsVarName, NumObs, Real2, &
                  num_levels = int(numchans, kind=integer64))
   do iChannel = 1, size(JediVarNamesWithChannels)
-    ! Retrieve data from JEDI
     call obsspace_get_db(ObsSpace, JediVarGroup, JediVarNamesWithChannels(iChannel), VarValue)
-
-    ! Fill the OPS data structures
-    if (localUseActualChans) then
-      where (VarValue /= MissingDouble)
-        Real2(:, Channels(iChannel)) = VarValue
-      end where
+    arrayindex = iChannel
+
+    ! if VAR channels have been assigned then jopa channels will be mapped to these var channels
+    ! Set up the size of the array, if channels are being pushed together an offset between
+    ! the var and jopa channel numbers is added onto the size of the array.
+    ! If not compressed the positions in the array are based on the actual channel number.
+
+    if (present(varChannels)) then
+      if (size(varChannels) > 0) then
+        if (iChannel <= size(varChannels)) then
+          if (compressChannels) then
+            offsetsize = abs(varChannels(1) - channels(1))
+            arrayindex = arrayindex + offsetsize
+          else
+            arrayindex = varChannels(iChannel)
+          end if
+        end if
+      else 
+        if (.not. compressChannels) then
+          arrayindex = Channels(iChannel)
+        end if
+      end if
     else
-     ! Fill the OPS data structures
-      where (VarValue /= MissingDouble)
-        Real2(:, iChannel+offset) = VarValue
-      end where
+      if (present(sizeOfVarobsArray)) then
+        if (sizeOfVarobsArray > size(channels)) then
+          if (.not. compressChannels) then
+            arrayindex = Channels(iChannel)
+          end if
+        end if
+      end if ! the end
     end if
+    where (VarValue /= MissingDouble)
+      Real2(:, arrayindex) = VarValue
+    end where
   end do
 end if ! Data not present? OPS will produce a warning -- we don't need to duplicate it.
+
 end subroutine opsinputs_fill_fillreal2d_norecords
 
 ! ------------------------------------------------------------------------------
@@ -1223,8 +1257,8 @@ end subroutine opsinputs_fill_fillreal2d_records
 !> We rely on warnings printed by the OPS code whenever data needed to output a requested varfield
 !> are not found.
 subroutine opsinputs_fill_fillreal2d( &
-  Hdr, OpsVarName, JediToOpsLayoutMapping, Real2, ObsSpace, Channels, VarobsLength, JediVarName, JediVarGroup, OffsetChans, &
-  useActualChans)
+  Hdr, OpsVarName, JediToOpsLayoutMapping, Real2, ObsSpace, Channels, &
+  VarobsLength, JediVarName, JediVarGroup, compressVarChannels, sizeOfVarobsArray, varChannels)
 implicit none
 
 ! Subroutine arguments:
@@ -1237,28 +1271,44 @@ subroutine opsinputs_fill_fillreal2d( &
 integer(integer64), intent(in)                      :: VarobsLength
 character(len=*), intent(in)                        :: JediVarName
 character(len=*), intent(in)                        :: JediVarGroup
-type(opsinputs_channeloffset), optional, intent(in) :: OffsetChans
-logical, optional, intent(in)                       :: useActualChans
+logical, optional, intent(in)                       :: compressVarChannels
+integer(integer64), optional, intent(in)            :: sizeOfVarobsArray
+integer(c_int), optional, intent(in)                :: varChannels(:)
+
+! local variables
+logical               :: compressChannels
+integer(integer64)    :: sizeOfVarobsArray_local
+integer(c_int), allocatable :: localvarChannels(:)
 
 ! Body:
+
+compressChannels = .true.
+if (present(compressVarChannels)) then
+  compressChannels = compressVarChannels
+end if
+
+if (present(varChannels)) then
+  allocate(localvarChannels(size(varChannels)))
+  localvarChannels = varChannels
+end if
+
+sizeOfVarobsArray_local = 0
+if (present(sizeOfVarobsArray)) then
+  sizeOfVarobsArray_local = sizeOfVarobsArray
+end if
+
 if (JediToOpsLayoutMapping % ConvertRecordsToMultilevelObs) then
   call opsinputs_fill_fillreal2d_records( &
     Hdr, OpsVarName, JediToOpsLayoutMapping, Real2, ObsSpace, VarobsLength, JediVarName, JediVarGroup)
 else
-  if (Present(OffsetChans)) then
-    if (Present(useActualChans)) then
-      call opsinputs_fill_fillreal2d_norecords( &
-        Hdr, OpsVarName, JediToOpsLayoutMapping % NumOpsObs, Real2, ObsSpace, Channels, &
-        JediVarName, JediVarGroup, OffsetChans, useActualChans)
-    else
-      call opsinputs_fill_fillreal2d_norecords( &
-        Hdr, OpsVarName, JediToOpsLayoutMapping % NumOpsObs, Real2, ObsSpace, Channels, &
-        JediVarName, JediVarGroup, OffsetChans)
-    end if
+  if (allocated(localvarChannels)) then
+    call opsinputs_fill_fillreal2d_norecords( &
+      Hdr, OpsVarName, JediToOpsLayoutMapping % NumOpsObs, Real2, ObsSpace, Channels, &
+      JediVarName, JediVarGroup, compressChannels, sizeOfVarobsArray_local, varChannels)
   else
     call opsinputs_fill_fillreal2d_norecords( &
       Hdr, OpsVarName, JediToOpsLayoutMapping % NumOpsObs, Real2, ObsSpace, Channels, &
-      JediVarName, JediVarGroup)
+      JediVarName, JediVarGroup, compressChannels, sizeOfVarobsArray_local)
   end if
 end if
 
@@ -2386,6 +2436,7 @@ function opsinputs_fill_varnames_with_channels(VarName, Channels) result(VarName
     write (VarNames(ichan),'(A,"_",I0)') VarName, Channels(ichan)
   end do
 end if
+
 end function opsinputs_fill_varnames_with_channels