From 39ac76f9ce4afe6efb65205c3f467b4881c3fe69 Mon Sep 17 00:00:00 2001
From: Brown <2mx@ornl.gov>
Date: Sat, 2 Apr 2022 11:21:12 -0400
Subject: [PATCH] add phys param writing

without changing object: CovarianceData
---
 sammy/src/fin/mfin0.f90                       |  5 ++-
 sammy/src/fin/mfin5.f90                       | 18 +++++++++--
 sammy/src/fit/mfit3.f90                       |  3 +-
 sammy/src/io/Hdf5IO.cpp                       | 29 ++++++++++++-----
 sammy/src/io/Hdf5IO.h                         | 10 +++++-
 sammy/src/io/interface/cix/Hdf5IO.cpp2f.xml   |  2 ++
 .../src/io/interface/cpp/Hdf5IOInterface.cpp  |  6 ++--
 sammy/src/io/interface/cpp/Hdf5IOInterface.h  |  4 +--
 sammy/src/io/interface/fortran/Hdf5IO_I.f90   |  6 ++--
 sammy/src/io/interface/fortran/Hdf5IO_M.f90   |  9 ++++--
 sammy/src/io/tests/Hdf5IOTest.cpp             |  4 ++-
 sammy/src/salmon/CovToolBox.cpp               | 31 +++++++++++++++++++
 sammy/src/salmon/CovToolBox.h                 | 12 +++++++
 13 files changed, 114 insertions(+), 25 deletions(-)

diff --git a/sammy/src/fin/mfin0.f90 b/sammy/src/fin/mfin0.f90
index bab45d6e2..f004ccd5c 100644
--- a/sammy/src/fin/mfin0.f90
+++ b/sammy/src/fin/mfin0.f90
@@ -6,6 +6,7 @@ module fin
 !
 ! *** Purpose -- Convert to physical parameters, output results
 !
+      ! ---- Common block modules ----------------
       use over_common_m
       use oops_common_m
       use fixedi_m
@@ -14,6 +15,7 @@ module fin
       use samxxx_common_m
       use oopsch_common_m
       use fixedr_m
+      ! ---- C++ interfaces ----------------------
       use SammyRMatrixParameters_M
       use SammyResonanceInfo_M
       use EndfData_common_m
@@ -21,11 +23,12 @@ module fin
       use ReadCovarianceInfo_m
       use ResonanceCovariance_M
       use AllocateFunctions_m
+      use SammyLptPrinting_m
+      ! ---- traditional modules -----------------
       use fin1
       use fin2
       use fin3
       use fin5
-      use SammyLptPrinting_m
       IMPLICIT DOUBLE PRECISION (a-h,o-z)
       real(kind=8),allocatable::Ratio(:)
       real(kind=8),allocatable,dimension(:)::A_Idum
diff --git a/sammy/src/fin/mfin5.f90 b/sammy/src/fin/mfin5.f90
index 334b12a96..8bba221f6 100644
--- a/sammy/src/fin/mfin5.f90
+++ b/sammy/src/fin/mfin5.f90
@@ -449,6 +449,7 @@ module fin5
 !
 ! *** PURPOSE -- WRITE COVARIANCE (ETC) FILE
 !
+      ! ---- common block modules ----------------
       use fixedi_m
       use ifwrit_m
       use samxxx_common_m
@@ -457,6 +458,7 @@ module fin5
       use namfil_common_m
       use par_parameter_names_common_m
       use EndfData_common_m
+      ! ---- C++ interfaces ----------------------
       use SammySpinGroupInfo_M
       use SammyIsoInfo_M
       use SammyResonanceInfo_M
@@ -464,10 +466,13 @@ module fin5
       use ReadCovarianceInfo_m
       use ResonanceCovariance_M
       use SammyRExternalInfo_M
+      use Hdf5IO_m
+
       IMPLICIT DOUBLE PRECISION (a-h,o-z)
-!
 !
       type(ResonanceCovariance)::physCov
+      type(Hdf5IO)::h5
+
       DIMENSION Parbrd(*), Iflbrd(*), &
          Pardet(*), Ifldet(*), Igrdet(*), &
          Polar(2,*),     &
@@ -511,11 +516,18 @@ module fin5
       integer,allocatable,dimension(:)::Ifliso
       integer,allocatable,dimension(:,:,:)::Iflext
       real(kind=8),allocatable,dimension(:,:,:)::Parext
-      integer::nradInfo
+      integer::nradInfo,wantUparam
       
       CHARACTER*3 Versn
       DATA Versn/'M7 '/
-!
+
+! ------------------------------------------------
+! Write out the RPCM to HDF5 in physical param space
+      call h5%initialize()
+      wantUparam = 0
+      call h5%writeRPCM(covData,resParData,trim(Samcov)//".h5",wantUparam)
+      call h5%destroy()
+! ------------------------------------------------
 !
       nradInfo = radFitFlags%getNumRadInfo()
       N = MAX0 (Nvpall, Nres*Ntotc2, Nres*Ntriag, &
diff --git a/sammy/src/fit/mfit3.f90 b/sammy/src/fit/mfit3.f90
index e48720633..9487ca02f 100644
--- a/sammy/src/fit/mfit3.f90
+++ b/sammy/src/fit/mfit3.f90
@@ -254,6 +254,7 @@ module fit3_m
          Engurr(*), Aold(3,*), Iflold(3,*), Ktold(*)
       type(ResonanceCovariance)::uCov, invCov
       type(Hdf5IO)::h5
+      integer, parameter :: wantUparam = 1
 !
 !
       CALL Newopn (62, Samcov, 1)
@@ -314,7 +315,7 @@ module fit3_m
       ! Write simple res param covariance matrix
       ! to HDF5 file (u-param & u-covar)
       ! -----------------------------------------
-      call h5%writeRPCM(covData,trim(Samcov)//".h5")
+      call h5%writeRPCM(covData,resParData,trim(Samcov)//".h5",wantUparam)
 
       RETURN
       END
diff --git a/sammy/src/io/Hdf5IO.cpp b/sammy/src/io/Hdf5IO.cpp
index baac8ff06..b65a43dfb 100644
--- a/sammy/src/io/Hdf5IO.cpp
+++ b/sammy/src/io/Hdf5IO.cpp
@@ -1,5 +1,4 @@
 #include "Hdf5IO.h"
-#include <cmath>
 
 namespace sammy {
     
@@ -102,23 +101,37 @@ namespace sammy {
         }
     }
 
-    void Hdf5IO::writeRPCM(CovarianceData & cov,std::string hdf5Filename){
+    void Hdf5IO::writeRPCM(CovarianceData & cov,
+                           SammyRMatrixParameters & resParData,
+                           std::string hdf5Filename,
+                           int rpcmIsReduced){
+        if( rpcmIsReduced < 0 || rpcmIsReduced > 1 ){
+            std::string errmsg("integer rpcmIsReduced must be 0 || 1");
+            throw std::runtime_error(errmsg);
+        }
         std::cout << "--------------------------------------------" << std::endl;
         std::cout << "writing resonance covariance to HDF5 file..." << std::endl;
         std::cout << "--------------------------------------------" << std::endl;
 
-        // get parameter covariance out
-        endf::ResonanceCovariance * parCov = cov.getUCovariance();
+        std::vector<double> param;
         int numParam = cov.getNumTotalParam();
+        endf::ResonanceCovariance * parCov;
+        // ---- get parameter covariance and parameters out ---------
+        if( rpcmIsReduced ){
+            parCov = cov.getUCovariance();
+            for( int i=0; i<numParam; ++i ) param.push_back(cov.getUParamValue(i));
+        }
+        else{
+            CovToolBox ctb;
+            param = ctb.getPhysParam(resParData);
+            parCov = cov.getCovariance();
+        }
+        // ---- check if length of params match covariance ----------
         int numRows = parCov->getNumRows();
         if( numParam != numRows ) {
             std::cerr << "numParam != numRows of covariance!" << std::endl;
         }
 
-        // get parameters out
-        std::vector<double> param;
-        for( int i=0; i<numParam; ++i ) param.push_back(cov.getUParamValue(i));
-
         double data[numRows][numRows];
         for( int i=0; i<numRows; ++i ){
             for( int j=0; j<numRows; ++j ){
diff --git a/sammy/src/io/Hdf5IO.h b/sammy/src/io/Hdf5IO.h
index 396439e51..fce378ccb 100644
--- a/sammy/src/io/Hdf5IO.h
+++ b/sammy/src/io/Hdf5IO.h
@@ -1,14 +1,19 @@
 #ifndef SAMMY_HDF5IO_H
 #define SAMMY_HDF5IO_H
 
+// ---- std ---------------------------
 #include <vector>
 #include <fstream>
 #include <iostream>
 #include <cstring>
 #include <stdexcept>
+#include <cmath>
+// ---- sammy -------------------------
 #include "OdfIO.h"
 #include "H5Cpp.h"
 #include "../endf/CovarianceData.h"
+#include "../endf/SammyRMatrixParameters.h"
+#include "../salmon/CovToolBox.h"
 
 namespace sammy{
 
@@ -60,7 +65,10 @@ namespace sammy{
          * @param cov class that contains all param. covariance info
          * @param hdf5Filename the string specifying the H5 file name
          ************************************************************/
-        void writeRPCM(CovarianceData & cov,std::string hdf5Filename);
+        void writeRPCM(CovarianceData & cov,
+                       SammyRMatrixParameters & resParData,
+                       std::string hdf5Filename,
+                       int rpcmIsReduced);
 
         /*************************************************************
          * Read the resonance parameter covariance matrix from an HDF5
diff --git a/sammy/src/io/interface/cix/Hdf5IO.cpp2f.xml b/sammy/src/io/interface/cix/Hdf5IO.cpp2f.xml
index 87b2cf43e..c86a09b1f 100644
--- a/sammy/src/io/interface/cix/Hdf5IO.cpp2f.xml
+++ b/sammy/src/io/interface/cix/Hdf5IO.cpp2f.xml
@@ -5,7 +5,9 @@
     <class name="Hdf5IO">
       <method name="writeRPCM">
          <param name="cov" type="CovarianceData"/>
+         <param name="resParData" type="SammyRMatrixParameters"/>
          <param name="hdf5Filename" type="char*" ctype="std::string"/>
+         <param name="rpcmIsReduced" type="int"/>
       </method>
       <method name="readRPCM">
          <param name="cov" type="CovarianceData"/>
diff --git a/sammy/src/io/interface/cpp/Hdf5IOInterface.cpp b/sammy/src/io/interface/cpp/Hdf5IOInterface.cpp
index 419d9310b..83d6fb957 100644
--- a/sammy/src/io/interface/cpp/Hdf5IOInterface.cpp
+++ b/sammy/src/io/interface/cpp/Hdf5IOInterface.cpp
@@ -2,16 +2,16 @@
 * This file has been dynamically generated by Class Interface Xml (CIX) 
 * DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 * If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Wed Jan 26 13:57:17 EST 2022
+* Date Generated: Sat Apr 02 11:08:38 EDT 2022
 * If any issues are experiences with this generated file that cannot be fixed
 * with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 */
 #include <string.h>
 #include "Hdf5IOInterface.h"
 using namespace sammy;
-void Hdf5IO_writeRPCM(void * Hdf5IO_ptr,CovarianceData * cov,char* hdf5Filename)
+void Hdf5IO_writeRPCM(void * Hdf5IO_ptr,CovarianceData * cov,SammyRMatrixParameters * resParData,char* hdf5Filename,int * rpcmIsReduced)
 {
-    ((Hdf5IO*)Hdf5IO_ptr)->writeRPCM(*cov,std::string(hdf5Filename));
+    ((Hdf5IO*)Hdf5IO_ptr)->writeRPCM(*cov,*resParData,std::string(hdf5Filename),*rpcmIsReduced);
 }
 
 void Hdf5IO_readRPCM(void * Hdf5IO_ptr,CovarianceData * cov,char* hdf5Filename)
diff --git a/sammy/src/io/interface/cpp/Hdf5IOInterface.h b/sammy/src/io/interface/cpp/Hdf5IOInterface.h
index b2ea36368..5e436ef2d 100644
--- a/sammy/src/io/interface/cpp/Hdf5IOInterface.h
+++ b/sammy/src/io/interface/cpp/Hdf5IOInterface.h
@@ -2,7 +2,7 @@
 * This file has been dynamically generated by Class Interface Xml (CIX) 
 * DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 * If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Wed Jan 26 13:57:17 EST 2022
+* Date Generated: Sat Apr 02 11:08:38 EDT 2022
 * If any issues are experiences with this generated file that cannot be fixed
 * with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 */
@@ -13,7 +13,7 @@ using namespace sammy;
 #ifdef __cplusplus
 extern "C" {
 #endif
-void Hdf5IO_writeRPCM(void * Hdf5IO_ptr,CovarianceData * cov,char* hdf5Filename);
+void Hdf5IO_writeRPCM(void * Hdf5IO_ptr,CovarianceData * cov,SammyRMatrixParameters * resParData,char* hdf5Filename,int * rpcmIsReduced);
 void Hdf5IO_readRPCM(void * Hdf5IO_ptr,CovarianceData * cov,char* hdf5Filename);
 void* Hdf5IO_initialize();
 void Hdf5IO_destroy(void * Hdf5IO_ptr);
diff --git a/sammy/src/io/interface/fortran/Hdf5IO_I.f90 b/sammy/src/io/interface/fortran/Hdf5IO_I.f90
index f9f92f7f4..ace498dff 100644
--- a/sammy/src/io/interface/fortran/Hdf5IO_I.f90
+++ b/sammy/src/io/interface/fortran/Hdf5IO_I.f90
@@ -2,19 +2,21 @@
 !! This file has been dynamically generated by Class Interface Xml (CIX) 
 !! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 !! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Wed Jan 26 13:57:17 EST 2022
+!! Date Generated: Sat Apr 02 11:08:38 EDT 2022
 !! If any issues are experiences with this generated file that cannot be fixed
 !! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 !!/
 module Hdf5IO_I
 use, intrinsic :: ISO_C_BINDING
 interface
-subroutine f_Hdf5IO_writeRPCM(Hdf5IO_ptr, cov,hdf5Filename ) BIND(C,name="Hdf5IO_writeRPCM")
+subroutine f_Hdf5IO_writeRPCM(Hdf5IO_ptr, cov,resParData,hdf5Filename,rpcmIsReduced ) BIND(C,name="Hdf5IO_writeRPCM")
     use,intrinsic :: ISO_C_BINDING
     implicit none
     type(C_PTR), value :: Hdf5IO_ptr;
     type(C_PTR), value :: cov;
+    type(C_PTR), value :: resParData;
     character(C_CHAR) :: hdf5Filename;
+    integer(C_INT) :: rpcmIsReduced;
 end subroutine
 subroutine f_Hdf5IO_readRPCM(Hdf5IO_ptr, cov,hdf5Filename ) BIND(C,name="Hdf5IO_readRPCM")
     use,intrinsic :: ISO_C_BINDING
diff --git a/sammy/src/io/interface/fortran/Hdf5IO_M.f90 b/sammy/src/io/interface/fortran/Hdf5IO_M.f90
index c6c492573..ad4efc070 100644
--- a/sammy/src/io/interface/fortran/Hdf5IO_M.f90
+++ b/sammy/src/io/interface/fortran/Hdf5IO_M.f90
@@ -2,7 +2,7 @@
 !! This file has been dynamically generated by Class Interface Xml (CIX) 
 !! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 !! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Wed Jan 26 13:57:17 EST 2022
+!! Date Generated: Sat Apr 02 11:08:38 EDT 2022
 !! If any issues are experiences with this generated file that cannot be fixed
 !! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 !!/
@@ -10,6 +10,7 @@ module Hdf5IO_M
 use, intrinsic :: ISO_C_BINDING
 use Hdf5IO_I
 use CovarianceData_M
+use SammyRMatrixParameters_M
 type Hdf5IO
     type(C_PTR) :: instance_ptr=C_NULL_PTR
     contains
@@ -19,12 +20,14 @@ type Hdf5IO
     procedure, pass(this) :: destroy => Hdf5IO_destroy
 end type Hdf5IO
 contains
-subroutine Hdf5IO_writeRPCM(this, cov, hdf5Filename)
+subroutine Hdf5IO_writeRPCM(this, cov, resParData, hdf5Filename, rpcmIsReduced)
     implicit none
     class(Hdf5IO)::this
     class(CovarianceData)::cov
+    class(SammyRMatrixParameters)::resParData
     character(len=*)::hdf5Filename
-    call f_Hdf5IO_writeRPCM(this%instance_ptr, cov%instance_ptr,hdf5Filename//C_NULL_CHAR)
+    integer(C_INT)::rpcmIsReduced
+    call f_Hdf5IO_writeRPCM(this%instance_ptr, cov%instance_ptr,resParData%instance_ptr,hdf5Filename//C_NULL_CHAR,rpcmIsReduced)
 end subroutine
 subroutine Hdf5IO_readRPCM(this, cov, hdf5Filename)
     implicit none
diff --git a/sammy/src/io/tests/Hdf5IOTest.cpp b/sammy/src/io/tests/Hdf5IOTest.cpp
index 53fb088ac..4a2b85100 100644
--- a/sammy/src/io/tests/Hdf5IOTest.cpp
+++ b/sammy/src/io/tests/Hdf5IOTest.cpp
@@ -9,8 +9,10 @@ TEST(Hdf5IO,readWriteRPCM)
 {
     sammy::CovarianceData covData;
     sammy::CovarianceData covDataR;
+    sammy::SammyRMatrixParameters resParData; //dummy
     sammy::Hdf5IO h5;
     sammy::Hdf5IO h5R;
+    int wantUparam = 1;
 
     // ----- set up the resonance covariance ---------------
     covData.makeNewUCovariance();
@@ -19,7 +21,7 @@ TEST(Hdf5IO,readWriteRPCM)
     cov->setCovariance(0,0,2.0); cov->setCovariance(0,1,1.0);
     cov->setCovariance(1,0,1.0); cov->setCovariance(1,1,4.0);
 
-    h5.writeRPCM(covData, "testreadwrite.h5");
+    h5.writeRPCM(covData,resParData, "testreadwrite.h5",wantUparam);
     h5R.readRPCM(covDataR,"testreadwrite.h5");
 
     endf::ResonanceCovariance * covR = covDataR.getUCovariance();
diff --git a/sammy/src/salmon/CovToolBox.cpp b/sammy/src/salmon/CovToolBox.cpp
index 9fb033497..0b2c7284e 100644
--- a/sammy/src/salmon/CovToolBox.cpp
+++ b/sammy/src/salmon/CovToolBox.cpp
@@ -48,4 +48,35 @@ namespace sammy{
             }
         }
     }
+
+    std::vector<double> CovToolBox::getPhysParam(SammyRMatrixParameters & resParData){
+
+        // --------------- grab all the *varied* pars --------------------------
+        int nr = resParData.getNumResonances();
+        std::vector<double> physParam;
+        for(int iRes =0;iRes<nr;++iRes){ // loop over resonances
+            sammy::SammyResonanceInfo *resInfo = resParData.getResonanceInfo(iRes);
+            int iGrp = resInfo->getSpinGroupIndex();
+            endf::RMatResonance *resonance = resParData.getResonance(resInfo);
+            sammy::SammySpinGroupInfo *spinInfo = resParData.getSpinGroupInfo(iGrp);
+
+            int parPerRes = spinInfo->getNumResPar();
+            for(int par=0;par<parPerRes;++par){ // loop over par's for each res
+                int covInd = 0; // >0 means varied or PUP'd
+                if(par==0){
+                    covInd = resInfo->getEnergyFitOption();
+                    if( covInd>0 && covInd!=3 ){
+                        physParam.push_back(resonance->getEres());
+                    }
+                }
+                else{
+                    covInd = resInfo->getChannelFitOption(par-1);
+                    if( covInd>0 && covInd!=3 ){
+                        physParam.push_back(resonance->getWidth(par-1));
+                    }
+                }
+            }
+        }
+        return physParam;
+    }
 }
\ No newline at end of file
diff --git a/sammy/src/salmon/CovToolBox.h b/sammy/src/salmon/CovToolBox.h
index 7a17de76f..b6a9df3a8 100644
--- a/sammy/src/salmon/CovToolBox.h
+++ b/sammy/src/salmon/CovToolBox.h
@@ -5,6 +5,9 @@
 
 #include "ScaleData/Core/io/CoverxFile.h"
 #include "ScaleUtils/EndfLib/CrossSectionCovariance.h"
+#include "endf/SammyRMatrixParameters.h"
+#include "endf/SammyResonanceInfo.h"
+#include "endf/CovarianceData.h"
 #include <iostream>
 #include <algorithm>
 
@@ -29,6 +32,15 @@ namespace sammy {
         void fillCrossCovFromCoverx(coverx::CoverxFile & cf, 
                                     endf::CrossCovList & cc);
 
+        /***********************************************************************
+         * Find varied parameters fill a vector with them and return the 
+         * physical parameters in the rigid "SAMMY order" in a vector of doubles
+         * 
+         * @param resParData the full r-matrix parameters object
+         * @return physParam vector of varied physical parameters
+         **********************************************************************/
+        std::vector<double> getPhysParam(SammyRMatrixParameters & resParData);
+
     private:
 
     };
-- 
GitLab