From 1deabfc9471fa344a3d0020f14d1f3dbf45e07c8 Mon Sep 17 00:00:00 2001
From: James E McClure <mcclurej@vt.edu>
Date: Wed, 6 Jun 2018 19:00:42 -0400
Subject: [PATCH] factored out uCT parameters

---
 analysis/uCT.cpp      | 13 +++++++------
 analysis/uCT.h        |  3 ++-
 tests/lbpm_uCT_pp.cpp |  8 ++++++--
 3 files changed, 15 insertions(+), 9 deletions(-)

diff --git a/analysis/uCT.cpp b/analysis/uCT.cpp
index 172de263..53c565b4 100644
--- a/analysis/uCT.cpp
+++ b/analysis/uCT.cpp
@@ -164,7 +164,8 @@ void solve( const Array<float>& VOL, Array<float>& Mean, Array<char>& ID,
 void refine( const Array<float>& Dist_coarse, 
     const Array<float>& VOL, Array<float>& Mean, Array<char>& ID,
     Array<float>& Dist, Array<float>& MultiScaleSmooth, Array<float>& NonLocalMean, 
-    fillHalo<float>& fillFloat, const Domain& Dm, int nprocx, int level )
+    fillHalo<float>& fillFloat, const Domain& Dm, int nprocx, int level,
+    float threshold, float lamda, float sigsq, int depth)
 {
     PROFILE_SCOPED(timer,"refine");
     int ratio[3] = { int(Dist.size(0)/Dist_coarse.size(0)),
@@ -175,7 +176,7 @@ void refine( const Array<float>& Dist_coarse,
     // Compute the median filter on the array and segment
     Med3D( VOL, Mean );
     fillFloat.fill( Mean );
-    segment( Mean, ID, 0.01 );
+    segment( Mean, ID, threshold );
     // If the ID has the wrong distance, set the distance to 0 and run a simple filter to set neighbors to 0
     for (size_t i=0; i<ID.length(); i++) {
         char id = Dist(i)>0 ? 1:0;
@@ -193,11 +194,11 @@ void refine( const Array<float>& Dist_coarse,
     Dist = imfilter::imfilter_separable<float>( Dist, {1,1,1}, filter_set, BC );
     fillFloat.fill( Dist );
     // Smooth the volume data
-    float lambda = 2*sqrt(double(ratio[0]*ratio[0]+ratio[1]*ratio[1]+ratio[2]*ratio[2]));
-    smooth( VOL, Dist, lambda, MultiScaleSmooth, fillFloat );
+    float h = 2*lamda*sqrt(double(ratio[0]*ratio[0]+ratio[1]*ratio[1]+ratio[2]*ratio[2]));
+    smooth( VOL, Dist, h, MultiScaleSmooth, fillFloat );
     // Compute non-local mean
-	int depth = 3;
-	float sigsq = 0.1;
+//	int depth = 3;
+//	float sigsq = 0.1;
 	int nlm_count = NLM3D( MultiScaleSmooth, Mean, Dist, NonLocalMean, depth, sigsq);
 	fillFloat.fill(NonLocalMean);
     segment( NonLocalMean, ID, 0.001 );
diff --git a/analysis/uCT.h b/analysis/uCT.h
index 32954d05..e1b1deca 100644
--- a/analysis/uCT.h
+++ b/analysis/uCT.h
@@ -31,7 +31,8 @@ void removeDisconnected( Array<char>& ID, const Domain& Dm );
 // Solve a level (without any coarse level information)
 void solve( const Array<float>& VOL, Array<float>& Mean, Array<char>& ID,
     Array<float>& Dist, Array<float>& MultiScaleSmooth, Array<float>& NonLocalMean, 
-    fillHalo<float>& fillFloat, const Domain& Dm, int nprocx );
+    fillHalo<float>& fillFloat, const Domain& Dm, int nprocx,
+    float threshold, float lamda, float sigsq, int depth);
 
 
 // Refine a solution from a coarse grid to a fine grid
diff --git a/tests/lbpm_uCT_pp.cpp b/tests/lbpm_uCT_pp.cpp
index a68233ce..1681e5aa 100644
--- a/tests/lbpm_uCT_pp.cpp
+++ b/tests/lbpm_uCT_pp.cpp
@@ -71,6 +71,10 @@ int main(int argc, char **argv)
 	int nprocz = nproc[2];
 
 	auto InputFile=uct_db->getScalar<std::string>( "InputFile" );    
+	auto rough_cutoff=uct_db->getScalar<float>( "rought_cutoff" );    
+	auto lamda=uct_db->getScalar<float>( "lamda" );    
+	auto nlm_sigsq=uct_db->getScalar<float>( "nlm_sigsq" );    
+	auto nlm_depth=uct_db->getScalar<int>( "nlm_depth" );    
 
 	//.......................................................................
 	// Reading the domain information file
@@ -243,7 +247,6 @@ int main(int argc, char **argv)
 		}
 	}
 
-
 	// Fill the source data for the coarse meshes
 	PROFILE_START("CoarsenMesh");
 	for (int i=1; i<N_levels; i++) {
@@ -269,7 +272,8 @@ int main(int argc, char **argv)
 	if (rank==0)
 		printf("Initialize full mesh\n");
 	solve( LOCVOL[0], Mean[0], ID[0], Dist[0], MultiScaleSmooth[0],
-			NonLocalMean[0], *fillFloat[0], *Dm[0], nprocx );
+			NonLocalMean[0], *fillFloat[0], *Dm[0], nprocx, 
+			rough_cutoff, lamda, nlm_sigsq, nlm_depth);
 	PROFILE_STOP("Solve fill mesh");
 	MPI_Barrier(comm);