Index: /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 26449)
+++ /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 26450)
@@ -125,4 +125,5 @@
 		case CalvingCrevasseDepthEnum:
 			parameters->AddObject(iomodel->CopyConstantObject("md.calving.crevasse_opening_stress",CalvingCrevasseDepthEnum));
+			parameters->AddObject(iomodel->CopyConstantObject("md.calving.crevasse_threshold",CalvingCrevasseThresholdEnum));
 			break;
 		case CalvingDev2Enum:
@@ -472,17 +473,12 @@
 
 	/*Intermediaries*/
-	int         calvinglaw;
-	IssmDouble  min_thickness,thickness,hab_fraction;
-	IssmDouble	crevassedepth,surface_crevasse,surface,critical_fraction;
-	IssmDouble  rho_ice,rho_water;
-	IssmDouble  bed,water_depth;
-	IssmDouble  levelset,sealevel;
-
+	int  calvinglaw;
 	femmodel->parameters->FindParam(&calvinglaw,CalvingLawEnum);
 
 	if(calvinglaw==CalvingMinthicknessEnum || calvinglaw==CalvingVonmisesEnum){
 
-		/*Get minimum thickness threshold*/
-		femmodel->parameters->FindParam(&min_thickness,CalvingMinthicknessEnum);
+		/*Intermediaries*/
+		IssmDouble thickness,bed,sealevel;
+		IssmDouble min_thickness = femmodel->parameters->FindParam(CalvingMinthicknessEnum);
 
 		/*Loop over all elements of this partition*/
@@ -518,4 +514,7 @@
    else if(calvinglaw==CalvingHabEnum){
 
+		/*Intermediaries*/
+		IssmDouble  thickness,water_depth,levelset,hab_fraction;
+
 		/*Get the fraction of the flotation thickness at the terminus*/
 		InputDuplicatex(femmodel,MaskIceLevelsetEnum,DistanceToCalvingfrontEnum);
@@ -526,6 +525,6 @@
 			Element* element  = xDynamicCast<Element*>(object);
 
-			rho_ice = element->FindParam(MaterialsRhoIceEnum);
-			rho_water = element->FindParam(MaterialsRhoSeawaterEnum);
+			IssmDouble rho_ice   = element->FindParam(MaterialsRhoIceEnum);
+			IssmDouble rho_water = element->FindParam(MaterialsRhoSeawaterEnum);
 
 			int      numnodes           = element->GetNumberOfNodes();
@@ -560,5 +559,7 @@
    else if(calvinglaw==CalvingCrevasseDepthEnum){
 
-		int   nflipped,local_nflipped;
+		/*Intermediaries*/
+		IssmDouble  levelset,crevassedepth,bed,surface_crevasse,thickness,surface;
+		int         nflipped,local_nflipped;
 		IssmDouble* constraint_nodes = NULL;
 
@@ -571,4 +572,6 @@
       int localmasters  = femmodel->nodes->NumberOfNodesLocal();
       Vector<IssmDouble>* vec_constraint_nodes = vec_constraint_nodes=new Vector<IssmDouble>(localmasters,numnodes);
+
+      IssmDouble crevasse_threshold = femmodel->parameters->FindParam(CalvingCrevasseThresholdEnum);
 
 		for(Object* & object : femmodel->elements->objects){
@@ -596,5 +599,5 @@
 					surface_input->GetInputValue(&surface,gauss);
 
-					if((surface_crevasse-surface>0. || crevassedepth-thickness>0.) && bed<0.){
+					if((surface_crevasse>surface || crevassedepth>crevasse_threshold*thickness) && bed<0.){
 						vec_constraint_nodes->SetValue(node->Pid(),1.0,INS_VAL);
 					}
@@ -616,10 +619,10 @@
 				Gauss*   gauss    = element->NewGauss();
 
-				Input*   levelset_input         = element->GetInput(DistanceToCalvingfrontEnum); _assert_(levelset_input);
-				Input*   crevassedepth_input    = element->GetInput(CrevasseDepthEnum); _assert_(crevassedepth_input);
-				Input*   bed_input              = element->GetInput(BedEnum); _assert_(bed_input);
-				Input*   surface_crevasse_input = element->GetInput(SurfaceCrevasseEnum); _assert_(surface_crevasse_input);
-				Input*   thickness_input        = element->GetInput(ThicknessEnum); _assert_(thickness_input);
-				Input*   surface_input          = element->GetInput(SurfaceEnum); _assert_(surface_input);
+				Input *levelset_input         = element->GetInput(DistanceToCalvingfrontEnum); _assert_(levelset_input);
+				Input *crevassedepth_input    = element->GetInput(CrevasseDepthEnum);          _assert_(crevassedepth_input);
+				Input *bed_input              = element->GetInput(BedEnum);                    _assert_(bed_input);
+				Input *surface_crevasse_input = element->GetInput(SurfaceCrevasseEnum);        _assert_(surface_crevasse_input);
+				Input *thickness_input        = element->GetInput(ThicknessEnum);              _assert_(thickness_input);
+				Input *surface_input          = element->GetInput(SurfaceEnum);                _assert_(surface_input);
 
 				/*Is this element connected to a node that should be calved*/
@@ -645,5 +648,5 @@
 						surface_input->GetInputValue(&surface,gauss);
 
-						if((surface_crevasse-surface>0. || crevassedepth-thickness>0.) && bed<0. && levelset>-300. && levelset<0. && constraint_nodes[node->Lid()]==0.){
+                  if((surface_crevasse>surface || crevassedepth>crevasse_threshold*thickness) && bed<0. && levelset>-300. && levelset<0. && constraint_nodes[node->Lid()]==0.){
 							local_nflipped++;
 							vec_constraint_nodes->SetValue(node->Pid(),1.0,INS_VAL);
@@ -655,5 +658,5 @@
 			/*Count how many new nodes were found*/
 			ISSM_MPI_Allreduce(&local_nflipped,&nflipped,1,ISSM_MPI_INT,ISSM_MPI_SUM,IssmComm::GetComm());
-			//_printf0_("Found "<<nflipped<<" to flip\n");
+			_printf0_("Found "<<nflipped<<" to flip\n");
 
 			/*Assemble and serialize flag vector*/
@@ -662,4 +665,5 @@
          femmodel->GetLocalVectorWithClonesNodes(&constraint_nodes,vec_constraint_nodes);
 		}
+
 		/*Free ressources:*/
 		delete vec_constraint_nodes;