Index: /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.cpp
===================================================================
--- /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.cpp	(revision 8607)
@@ -11,21 +11,33 @@
 #include "../Responsex/Responsex.h"
 
-void CostFunctionx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,int response){
+void CostFunctionx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters){
 
 	/*Intermediary*/
 	int    i;
-	double fit;
+	int     num_responses;
 	double S;
 	Element* element=NULL;
+	int*     responses=NULL;
 
 	/*output: */
-	double J;
+	double J,Jplus;
+
+	
+	/*Recover parameters*/
+	parameters->FindParam(&num_responses,NumResponsesEnum);
+	parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
+
+	/*Get response*/
+	J=0;
+	for(int i=0;i<num_responses;i++){
+		Responsex(&Jplus,elements,nodes,vertices,loads,materials,parameters,EnumToStringx(responses[i]),false,i); //False means DO NOT process units
+		J+=Jplus;
+	}
+
+	/*REST TO BE DELETED*/
+
+	/*Add Regularization terms: */
 	double Jreg=0;
 	double Jreg_sum;
-	
-	/*Get response*/
-	Responsex(&J,elements,nodes,vertices,loads,materials,parameters,EnumToStringx(response),false); //False means DO NOT process units
-
-	/*Add Regularization terms: */
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
@@ -39,4 +51,5 @@
 
 	/*Assign output pointers: */
+	xfree((void**)&responses);
 	*pJ=J;
 }
Index: /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.h
===================================================================
--- /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/CostFunctionx/CostFunctionx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void CostFunctionx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,int response);
+void CostFunctionx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters);
 
 #endif
Index: /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp
===================================================================
--- /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp	(revision 8607)
@@ -76,5 +76,5 @@
 			
 			//Responsex(responses_pointer,elements,nodes, vertices,loads,materials, parameters,root,process_units);
-			Responsex(&femmodel_response,elements,nodes, vertices,loads,materials, parameters,root,process_units);
+			Responsex(&femmodel_response,elements,nodes, vertices,loads,materials, parameters,root,process_units,0);//0 is the index for weights
 			
 			if(my_rank==0){
@@ -95,5 +95,5 @@
 			
 			/*perfectly normal response function: */
-			Responsex(&femmodel_response,elements,nodes, vertices,loads,materials, parameters,root,process_units);
+			Responsex(&femmodel_response,elements,nodes, vertices,loads,materials, parameters,root,process_units,0);//0 is the weight index
 
 			if(my_rank==0){
Index: /issm/trunk/src/c/modules/Responsex/Responsex.cpp
===================================================================
--- /issm/trunk/src/c/modules/Responsex/Responsex.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/Responsex/Responsex.cpp	(revision 8607)
@@ -16,26 +16,26 @@
 #include "../modules.h"
 
-void Responsex(double* responses,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,char* response_descriptor,bool process_units){
+void Responsex(double* responses,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,char* response_descriptor,bool process_units,int weight_index){
 
 	switch (StringToEnumx(response_descriptor)){
 		
-		case MinVelEnum: MinVelx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case MaxVelEnum: MaxVelx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case MinVxEnum:  MinVxx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);  break;
-		case MaxVxEnum:  MaxVxx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);  break;
-		case MaxAbsVxEnum: MaxAbsVxx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case MinVyEnum:  MinVyx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);  break;
-		case MaxVyEnum:  MaxVyx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);  break;
-		case MaxAbsVyEnum: MaxAbsVyx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case MinVzEnum: MinVzx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);   break;
-		case MaxVzEnum: MaxVzx( responses, elements,nodes, vertices, loads, materials, parameters,process_units);   break;
-		case MaxAbsVzEnum: MaxAbsVzx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case MassFluxEnum: MassFluxx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case SurfaceAbsVelMisfitEnum: SurfaceAbsVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case SurfaceRelVelMisfitEnum: SurfaceRelVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case SurfaceLogVelMisfitEnum: SurfaceLogVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case SurfaceLogVxVyMisfitEnum: SurfaceLogVxVyMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case SurfaceAverageVelMisfitEnum: SurfaceAverageVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
-		case ThicknessAbsMisfitEnum: ThicknessAbsMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MinVelEnum:                 MinVelx(                  responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxVelEnum:                 MaxVelx(                  responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MinVxEnum:                  MinVxx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxVxEnum:                  MaxVxx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxAbsVxEnum:               MaxAbsVxx(                responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MinVyEnum:                  MinVyx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxVyEnum:                  MaxVyx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxAbsVyEnum:               MaxAbsVyx(                responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MinVzEnum:                  MinVzx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxVzEnum:                  MaxVzx(                   responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MaxAbsVzEnum:               MaxAbsVzx(                responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case MassFluxEnum:               MassFluxx(                responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
+		case SurfaceAbsVelMisfitEnum:    SurfaceAbsVelMisfitx(     responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
+		case SurfaceRelVelMisfitEnum:    SurfaceRelVelMisfitx(     responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
+		case SurfaceLogVelMisfitEnum:    SurfaceLogVelMisfitx(     responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
+		case SurfaceLogVxVyMisfitEnum:   SurfaceLogVxVyMisfitx(    responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
+		case SurfaceAverageVelMisfitEnum:SurfaceAverageVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
+		case ThicknessAbsMisfitEnum:     ThicknessAbsMisfitx(      responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break;
 		default: _error_(" response descriptor \"%s\" not supported yet!",response_descriptor); break;
 	}
Index: /issm/trunk/src/c/modules/Responsex/Responsex.h
===================================================================
--- /issm/trunk/src/c/modules/Responsex/Responsex.h	(revision 8606)
+++ /issm/trunk/src/c/modules/Responsex/Responsex.h	(revision 8607)
@@ -9,5 +9,5 @@
 #include "../../Container/Container.h"
 
-void Responsex(double* presponse,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,char* response_descriptor,bool process_units);
+void Responsex(double* presponse,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,char* response_descriptor,bool process_units,int weight_index);
 
 #endif  /* _RESPONSESXX_H */
Index: /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.cpp	(revision 8607)
@@ -10,5 +10,5 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void SurfaceAbsVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void SurfaceAbsVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units,int weight_index){
 
 	/*Intermediary*/
@@ -23,5 +23,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->SurfaceAbsVelMisfit(process_units);
+		J+=element->SurfaceAbsVelMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void SurfaceAbsVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void SurfaceAbsVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.cpp	(revision 8607)
@@ -11,5 +11,5 @@
 #include "../SurfaceAreax/SurfaceAreax.h"
 
-void SurfaceAverageVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void SurfaceAverageVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units,int weight_index){
 
 	/*Intermediary*/
@@ -27,5 +27,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->SurfaceAverageVelMisfit(process_units);
+		J+=element->SurfaceAverageVelMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void SurfaceAverageVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void SurfaceAverageVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.cpp	(revision 8607)
@@ -10,5 +10,5 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void SurfaceLogVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void SurfaceLogVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units,int weight_index){
 
 	/*Intermediary*/
@@ -23,5 +23,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->SurfaceLogVelMisfit(process_units);
+		J+=element->SurfaceLogVelMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void SurfaceLogVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void SurfaceLogVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.cpp	(revision 8607)
@@ -10,5 +10,5 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void SurfaceLogVxVyMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void SurfaceLogVxVyMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units,int weight_index){
 
 	/*Intermediary*/
@@ -23,5 +23,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->SurfaceLogVxVyMisfit(process_units);
+		J+=element->SurfaceLogVxVyMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void SurfaceLogVxVyMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void SurfaceLogVxVyMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.cpp	(revision 8607)
@@ -10,5 +10,5 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void SurfaceRelVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void SurfaceRelVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units,int weight_index){
 
 	/*Intermediary*/
@@ -23,5 +23,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->SurfaceRelVelMisfit(process_units);
+		J+=element->SurfaceRelVelMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void SurfaceRelVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void SurfaceRelVelMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.cpp
===================================================================
--- /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.cpp	(revision 8606)
+++ /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.cpp	(revision 8607)
@@ -10,5 +10,5 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void ThicknessAbsMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units){
+void ThicknessAbsMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials,Parameters* parameters,bool process_units, int weight_index){
 
 	/*Intermediary*/
@@ -23,5 +23,5 @@
 	for (i=0;i<elements->Size();i++){
 		element=(Element*)elements->GetObjectByOffset(i);
-		J+=element->ThicknessAbsMisfit(process_units);
+		J+=element->ThicknessAbsMisfit(process_units,weight_index);
 	}
 
Index: /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.h
===================================================================
--- /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.h	(revision 8606)
+++ /issm/trunk/src/c/modules/ThicknessAbsMisfitx/ThicknessAbsMisfitx.h	(revision 8607)
@@ -10,5 +10,5 @@
 
 /* local prototypes: */
-void ThicknessAbsMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units);
+void ThicknessAbsMisfitx( double* pJ, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,bool process_units,int weight_index);
 
 #endif
Index: /issm/trunk/src/c/objects/Elements/Element.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Element.h	(revision 8606)
+++ /issm/trunk/src/c/objects/Elements/Element.h	(revision 8607)
@@ -43,10 +43,10 @@
 		virtual void   GradjDrag(Vec gradient)=0;
 		virtual void   GradjB(Vec gradient)=0;
-		virtual double ThicknessAbsMisfit(bool process_units)=0;
-		virtual double SurfaceAbsVelMisfit(bool process_units)=0;
-		virtual double SurfaceRelVelMisfit(bool process_units)=0;
-		virtual double SurfaceLogVelMisfit(bool process_units)=0;
-		virtual double SurfaceLogVxVyMisfit(bool process_units)=0;
-		virtual double SurfaceAverageVelMisfit(bool process_units)=0;
+		virtual double ThicknessAbsMisfit(bool process_units  ,int weight_index)=0;
+		virtual double SurfaceAbsVelMisfit(bool process_units ,int weight_index)=0;
+		virtual double SurfaceRelVelMisfit(bool process_units ,int weight_index)=0;
+		virtual double SurfaceLogVelMisfit(bool process_units ,int weight_index)=0;
+		virtual double SurfaceLogVxVyMisfit(bool process_units,int weight_index)=0;
+		virtual double SurfaceAverageVelMisfit(bool process_units,int weight_index)=0;
 		virtual double RegularizeInversion(void)=0;
 		virtual double SurfaceArea(void)=0;
Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 8606)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 8607)
@@ -6845,5 +6845,5 @@
 /*}}}*/
 /*FUNCTION Penta::SurfaceAverageVelMisfit {{{1*/
-double Penta::SurfaceAverageVelMisfit(bool process_units){
+double Penta::SurfaceAverageVelMisfit(bool process_units,int weight_index){
 
 	int    approximation;
@@ -6868,5 +6868,5 @@
 		 * and compute SurfaceAverageVelMisfit*/
 		tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
-		J=tria->SurfaceAverageVelMisfit(process_units);
+		J=tria->SurfaceAverageVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6875,5 +6875,5 @@
 
 		tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-		J=tria->SurfaceAverageVelMisfit(process_units);
+		J=tria->SurfaceAverageVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6882,5 +6882,5 @@
 /*}}}*/
 /*FUNCTION Penta::SurfaceAbsVelMisfit {{{1*/
-double Penta::SurfaceAbsVelMisfit(bool process_units){
+double Penta::SurfaceAbsVelMisfit(bool process_units,int weight_index){
 
 	int    approximation;
@@ -6905,5 +6905,5 @@
 		 * and compute SurfaceAbsVelMisfit*/
 		tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
-		J=tria->SurfaceAbsVelMisfit(process_units);
+		J=tria->SurfaceAbsVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6912,5 +6912,5 @@
 
 		tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-		J=tria->SurfaceAbsVelMisfit(process_units);
+		J=tria->SurfaceAbsVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6919,5 +6919,5 @@
 /*}}}*/
 /*FUNCTION Penta::SurfaceLogVelMisfit {{{1*/
-double Penta::SurfaceLogVelMisfit(bool process_units){
+double Penta::SurfaceLogVelMisfit(bool process_units,int weight_index){
 
 	int    approximation;
@@ -6942,5 +6942,5 @@
 		 * and compute SurfaceLogVelMisfit*/
 		tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
-		J=tria->SurfaceLogVelMisfit(process_units);
+		J=tria->SurfaceLogVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6949,5 +6949,5 @@
 
 		tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-		J=tria->SurfaceLogVelMisfit(process_units);
+		J=tria->SurfaceLogVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6956,5 +6956,5 @@
 /*}}}*/
 /*FUNCTION Penta::SurfaceLogVxVyMisfit {{{1*/
-double Penta::SurfaceLogVxVyMisfit(bool process_units){
+double Penta::SurfaceLogVxVyMisfit(bool process_units,int weight_index){
 
 	double J;
@@ -6981,5 +6981,5 @@
 		 * and compute SurfaceLogVxVyMisfit*/
 		tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
-		J=tria->SurfaceLogVxVyMisfit(process_units);
+		J=tria->SurfaceLogVxVyMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -6988,5 +6988,5 @@
 
 		tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-		J=tria->SurfaceLogVxVyMisfit(process_units);
+		J=tria->SurfaceLogVxVyMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -7019,5 +7019,5 @@
 /*}}}*/
 /*FUNCTION Penta::SurfaceRelVelMisfit {{{1*/
-double Penta::SurfaceRelVelMisfit(bool process_units){
+double Penta::SurfaceRelVelMisfit(bool process_units,int weight_index){
 
 	int    approximation;
@@ -7042,5 +7042,5 @@
 		 * and compute SurfaceRelVelMisfit*/
 		tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
-		J=tria->SurfaceRelVelMisfit(process_units);
+		J=tria->SurfaceRelVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -7049,5 +7049,5 @@
 
 		tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-		J=tria->SurfaceRelVelMisfit(process_units);
+		J=tria->SurfaceRelVelMisfit(process_units,weight_index);
 		delete tria->matice; delete tria;
 		return J;
@@ -7100,5 +7100,5 @@
 }
 /*FUNCTION Penta::ThicknessAbsMisfit {{{1*/
-double Penta::ThicknessAbsMisfit(bool process_units){
+double Penta::ThicknessAbsMisfit(bool process_units,int weight_index){
 
 	int    approximation;
@@ -7113,6 +7113,6 @@
 	_error_("Not implemented yet");
 
-	tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
-	J=tria->ThicknessAbsMisfit(process_units);
+	tria=(Tria*)SpawnTria(0,1,2);
+	J=tria->ThicknessAbsMisfit(process_units,weight_index);
 	delete tria->matice; delete tria;
 	return J;
Index: /issm/trunk/src/c/objects/Elements/Penta.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.h	(revision 8606)
+++ /issm/trunk/src/c/objects/Elements/Penta.h	(revision 8607)
@@ -118,10 +118,10 @@
 		void   MinVy(double* pminvy, bool process_units);
 		void   MinVz(double* pminvz, bool process_units);
-		double ThicknessAbsMisfit(bool process_units);
-		double SurfaceAbsVelMisfit(bool process_units);
-		double SurfaceRelVelMisfit(bool process_units);
-		double SurfaceLogVelMisfit(bool process_units);
-		double SurfaceLogVxVyMisfit(bool process_units);
-		double SurfaceAverageVelMisfit(bool process_units);
+		double ThicknessAbsMisfit(     bool process_units,int weight_index);
+		double SurfaceAbsVelMisfit(    bool process_units,int weight_index);
+		double SurfaceRelVelMisfit(    bool process_units,int weight_index);
+		double SurfaceLogVelMisfit(    bool process_units,int weight_index);
+		double SurfaceLogVxVyMisfit(   bool process_units,int weight_index);
+		double SurfaceAverageVelMisfit(bool process_units,int weight_index);
 		void   PatchFill(int* pcount, Patch* patch);
 		void   PatchSize(int* pnumrows, int* pnumvertices,int* pnumnodes);
Index: /issm/trunk/src/c/objects/Elements/Tria.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 8606)
+++ /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 8607)
@@ -1706,5 +1706,4 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -1714,103 +1713,107 @@
 
 		/*Loop over all requested responses*/
-		for(resp=0;resp<num_responses;resp++) switch(responses[resp]){
-
-			case SurfaceAbsVelMisfitEnum:
-				/*
-				 *      1  [           2              2 ]
-				 * J = --- | (u - u   )  +  (v - v   )  |
-				 *      2  [       obs            obs   ]
-				 *
-				 *        dJ
-				 * DU = - -- = (u   - u )
-				 *        du     obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					dux=vxobs-vx;
-					duy=vyobs-vy;
-					pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceRelVelMisfitEnum:
-				/*
-				 *      1  [     \bar{v}^2             2   \bar{v}^2              2 ]
-				 * J = --- | -------------  (u - u   ) + -------------  (v - v   )  |
-				 *      2  [  (u   + eps)^2       obs    (v   + eps)^2       obs    ]
-				 *              obs                        obs                      
-				 *
-				 *        dJ     \bar{v}^2
-				 * DU = - -- = ------------- (u   - u )
-				 *        du   (u   + eps)^2    obs
-				 *               obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					scalex=pow(meanvel/(vxobs+epsvel),2.); if(vxobs==0)scalex=0;
-					scaley=pow(meanvel/(vyobs+epsvel),2.); if(vyobs==0)scaley=0;
-					dux=scalex*(vxobs-vx);
-					duy=scaley*(vyobs-vy);
-					pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceLogVelMisfitEnum:
-				/*
-				 *                 [        vel + eps     ] 2
-				 * J = 4 \bar{v}^2 | log ( -----------  ) |  
-				 *                 [       vel   + eps    ]
-				 *                            obs
-				 *
-				 *        dJ                 2 * log(...)
-				 * DU = - -- = - 4 \bar{v}^2 -------------  u
-				 *        du                 vel^2 + eps
-				 *            
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					velocity_mag    =sqrt(pow(vx,   2.)+pow(vy,   2.))+epsvel;
-					obs_velocity_mag=sqrt(pow(vxobs,2.)+pow(vyobs,2.))+epsvel;
-					scale=-8*pow(meanvel,2.)/pow(velocity_mag,2.)*log(velocity_mag/obs_velocity_mag);
-					dux=scale*vx;
-					duy=scale*vy;
-					pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceAverageVelMisfitEnum:
-				/*
-				 *      1                    2              2
-				 * J = ---  sqrt(  (u - u   )  +  (v - v   )  )
-				 *      S                obs            obs
-				 *
-				 *        dJ      1       1 
-				 * DU = - -- = - --- ----------- * 2 (u - u   )
-				 *        du      S  2 sqrt(...)           obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					scale=1./(S*2*sqrt(pow(vx-vxobs,2.)+pow(vy-vyobs,2.))+epsvel);
-					dux=scale*(vxobs-vx);
-					duy=scale*(vyobs-vy);
-					pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceLogVxVyMisfitEnum:
-				/*
-				 *      1            [        |u| + eps     2          |v| + eps     2  ]
-				 * J = --- \bar{v}^2 | log ( -----------  )   +  log ( -----------  )   |  
-				 *      2            [       |u    |+ eps              |v    |+ eps     ]
-				 *                              obs                       obs
-				 *        dJ                              1      u                             1
-				 * DU = - -- = - \bar{v}^2 log(u...) --------- ----  ~ - \bar{v}^2 log(u...) ------
-				 *        du                         |u| + eps  |u|                           u + eps
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					dux = - pow(meanvel,2.) * log((fabs(vx)+epsvel)/(fabs(vxobs)+epsvel)) / (vx+epsvel);
-					duy = - pow(meanvel,2.) * log((fabs(vy)+epsvel)/(fabs(vyobs)+epsvel)) / (vy+epsvel);
-					pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			default:
-				_error_("response %s not supported yet",EnumToStringx(responses[resp]));
+		for(resp=0;resp<num_responses;resp++){
+
+			weights_input->GetParameterValue(&weight,gauss,resp);
+
+			switch(responses[resp]){
+				case SurfaceAbsVelMisfitEnum:
+					/*
+					 *      1  [           2              2 ]
+					 * J = --- | (u - u   )  +  (v - v   )  |
+					 *      2  [       obs            obs   ]
+					 *
+					 *        dJ
+					 * DU = - -- = (u   - u )
+					 *        du     obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						dux=vxobs-vx;
+						duy=vyobs-vy;
+						pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceRelVelMisfitEnum:
+					/*
+					 *      1  [     \bar{v}^2             2   \bar{v}^2              2 ]
+					 * J = --- | -------------  (u - u   ) + -------------  (v - v   )  |
+					 *      2  [  (u   + eps)^2       obs    (v   + eps)^2       obs    ]
+					 *              obs                        obs                      
+					 *
+					 *        dJ     \bar{v}^2
+					 * DU = - -- = ------------- (u   - u )
+					 *        du   (u   + eps)^2    obs
+					 *               obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						scalex=pow(meanvel/(vxobs+epsvel),2.); if(vxobs==0)scalex=0;
+						scaley=pow(meanvel/(vyobs+epsvel),2.); if(vyobs==0)scaley=0;
+						dux=scalex*(vxobs-vx);
+						duy=scaley*(vyobs-vy);
+						pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceLogVelMisfitEnum:
+					/*
+					 *                 [        vel + eps     ] 2
+					 * J = 4 \bar{v}^2 | log ( -----------  ) |  
+					 *                 [       vel   + eps    ]
+					 *                            obs
+					 *
+					 *        dJ                 2 * log(...)
+					 * DU = - -- = - 4 \bar{v}^2 -------------  u
+					 *        du                 vel^2 + eps
+					 *            
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						velocity_mag    =sqrt(pow(vx,   2.)+pow(vy,   2.))+epsvel;
+						obs_velocity_mag=sqrt(pow(vxobs,2.)+pow(vyobs,2.))+epsvel;
+						scale=-8*pow(meanvel,2.)/pow(velocity_mag,2.)*log(velocity_mag/obs_velocity_mag);
+						dux=scale*vx;
+						duy=scale*vy;
+						pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceAverageVelMisfitEnum:
+					/*
+					 *      1                    2              2
+					 * J = ---  sqrt(  (u - u   )  +  (v - v   )  )
+					 *      S                obs            obs
+					 *
+					 *        dJ      1       1 
+					 * DU = - -- = - --- ----------- * 2 (u - u   )
+					 *        du      S  2 sqrt(...)           obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						scale=1./(S*2*sqrt(pow(vx-vxobs,2.)+pow(vy-vyobs,2.))+epsvel);
+						dux=scale*(vxobs-vx);
+						duy=scale*(vyobs-vy);
+						pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceLogVxVyMisfitEnum:
+					/*
+					 *      1            [        |u| + eps     2          |v| + eps     2  ]
+					 * J = --- \bar{v}^2 | log ( -----------  )   +  log ( -----------  )   |  
+					 *      2            [       |u    |+ eps              |v    |+ eps     ]
+					 *                              obs                       obs
+					 *        dJ                              1      u                             1
+					 * DU = - -- = - \bar{v}^2 log(u...) --------- ----  ~ - \bar{v}^2 log(u...) ------
+					 *        du                         |u| + eps  |u|                           u + eps
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						dux = - pow(meanvel,2.) * log((fabs(vx)+epsvel)/(fabs(vxobs)+epsvel)) / (vx+epsvel);
+						duy = - pow(meanvel,2.) * log((fabs(vy)+epsvel)/(fabs(vyobs)+epsvel)) / (vy+epsvel);
+						pe->values[i*NDOF2+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF2+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				default:
+					_error_("response %s not supported yet",EnumToStringx(responses[resp]));
+			}
 		}
 	}
@@ -1870,5 +1873,4 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -1878,103 +1880,108 @@
 
 		/*Loop over all requested responses*/
-		for(resp=0;resp<num_responses;resp++) switch(responses[resp]){
-
-			case SurfaceAbsVelMisfitEnum:
-				/*
-				 *      1  [           2              2 ]
-				 * J = --- | (u - u   )  +  (v - v   )  |
-				 *      2  [       obs            obs   ]
-				 *
-				 *        dJ
-				 * DU = - -- = (u   - u )
-				 *        du     obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					dux=vxobs-vx;
-					duy=vyobs-vy;
-					pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceRelVelMisfitEnum:
-				/*
-				 *      1  [     \bar{v}^2             2   \bar{v}^2              2 ]
-				 * J = --- | -------------  (u - u   ) + -------------  (v - v   )  |
-				 *      2  [  (u   + eps)^2       obs    (v   + eps)^2       obs    ]
-				 *              obs                        obs                      
-				 *
-				 *        dJ     \bar{v}^2
-				 * DU = - -- = ------------- (u   - u )
-				 *        du   (u   + eps)^2    obs
-				 *               obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					scalex=pow(meanvel/(vxobs+epsvel),2.); if(vxobs==0)scalex=0;
-					scaley=pow(meanvel/(vyobs+epsvel),2.); if(vyobs==0)scaley=0;
-					dux=scalex*(vxobs-vx);
-					duy=scaley*(vyobs-vy);
-					pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceLogVelMisfitEnum:
-				/*
-				 *                 [        vel + eps     ] 2
-				 * J = 4 \bar{v}^2 | log ( -----------  ) |  
-				 *                 [       vel   + eps    ]
-				 *                            obs
-				 *
-				 *        dJ                 2 * log(...)
-				 * DU = - -- = - 4 \bar{v}^2 -------------  u
-				 *        du                 vel^2 + eps
-				 *            
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					velocity_mag    =sqrt(pow(vx,   2.)+pow(vy,   2.))+epsvel;
-					obs_velocity_mag=sqrt(pow(vxobs,2.)+pow(vyobs,2.))+epsvel;
-					scale=-8*pow(meanvel,2.)/pow(velocity_mag,2.)*log(velocity_mag/obs_velocity_mag);
-					dux=scale*vx;
-					duy=scale*vy;
-					pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceAverageVelMisfitEnum:
-				/*
-				 *      1                    2              2
-				 * J = ---  sqrt(  (u - u   )  +  (v - v   )  )
-				 *      S                obs            obs
-				 *
-				 *        dJ      1       1 
-				 * DU = - -- = - --- ----------- * 2 (u - u   )
-				 *        du      S  2 sqrt(...)           obs
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					scale=1./(S*2*sqrt(pow(vx-vxobs,2.)+pow(vy-vyobs,2.))+epsvel);
-					dux=scale*(vxobs-vx);
-					duy=scale*(vyobs-vy);
-					pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			case SurfaceLogVxVyMisfitEnum:
-				/*
-				 *      1            [        |u| + eps     2          |v| + eps     2  ]
-				 * J = --- \bar{v}^2 | log ( -----------  )   +  log ( -----------  )   |  
-				 *      2            [       |u    |+ eps              |v    |+ eps     ]
-				 *                              obs                       obs
-				 *        dJ                              1      u                             1
-				 * DU = - -- = - \bar{v}^2 log(u...) --------- ----  ~ - \bar{v}^2 log(u...) ------
-				 *        du                         |u| + eps  |u|                           u + eps
-				 */
-				for (i=0;i<NUMVERTICES;i++){
-					dux = - pow(meanvel,2.) * log((fabs(vx)+epsvel)/(fabs(vxobs)+epsvel)) / (vx+epsvel);
-					duy = - pow(meanvel,2.) * log((fabs(vy)+epsvel)/(fabs(vyobs)+epsvel)) / (vy+epsvel);
-					pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
-					pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
-				}
-				break;
-			default:
-				_error_("response %s not supported yet",EnumToStringx(responses[resp]));
+		for(resp=0;resp<num_responses;resp++){
+
+			weights_input->GetParameterValue(&weight,gauss,resp);
+
+			switch(responses[resp]){
+
+				case SurfaceAbsVelMisfitEnum:
+					/*
+					 *      1  [           2              2 ]
+					 * J = --- | (u - u   )  +  (v - v   )  |
+					 *      2  [       obs            obs   ]
+					 *
+					 *        dJ
+					 * DU = - -- = (u   - u )
+					 *        du     obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						dux=vxobs-vx;
+						duy=vyobs-vy;
+						pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceRelVelMisfitEnum:
+					/*
+					 *      1  [     \bar{v}^2             2   \bar{v}^2              2 ]
+					 * J = --- | -------------  (u - u   ) + -------------  (v - v   )  |
+					 *      2  [  (u   + eps)^2       obs    (v   + eps)^2       obs    ]
+					 *              obs                        obs                      
+					 *
+					 *        dJ     \bar{v}^2
+					 * DU = - -- = ------------- (u   - u )
+					 *        du   (u   + eps)^2    obs
+					 *               obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						scalex=pow(meanvel/(vxobs+epsvel),2.); if(vxobs==0)scalex=0;
+						scaley=pow(meanvel/(vyobs+epsvel),2.); if(vyobs==0)scaley=0;
+						dux=scalex*(vxobs-vx);
+						duy=scaley*(vyobs-vy);
+						pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceLogVelMisfitEnum:
+					/*
+					 *                 [        vel + eps     ] 2
+					 * J = 4 \bar{v}^2 | log ( -----------  ) |  
+					 *                 [       vel   + eps    ]
+					 *                            obs
+					 *
+					 *        dJ                 2 * log(...)
+					 * DU = - -- = - 4 \bar{v}^2 -------------  u
+					 *        du                 vel^2 + eps
+					 *            
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						velocity_mag    =sqrt(pow(vx,   2.)+pow(vy,   2.))+epsvel;
+						obs_velocity_mag=sqrt(pow(vxobs,2.)+pow(vyobs,2.))+epsvel;
+						scale=-8*pow(meanvel,2.)/pow(velocity_mag,2.)*log(velocity_mag/obs_velocity_mag);
+						dux=scale*vx;
+						duy=scale*vy;
+						pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceAverageVelMisfitEnum:
+					/*
+					 *      1                    2              2
+					 * J = ---  sqrt(  (u - u   )  +  (v - v   )  )
+					 *      S                obs            obs
+					 *
+					 *        dJ      1       1 
+					 * DU = - -- = - --- ----------- * 2 (u - u   )
+					 *        du      S  2 sqrt(...)           obs
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						scale=1./(S*2*sqrt(pow(vx-vxobs,2.)+pow(vy-vyobs,2.))+epsvel);
+						dux=scale*(vxobs-vx);
+						duy=scale*(vyobs-vy);
+						pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				case SurfaceLogVxVyMisfitEnum:
+					/*
+					 *      1            [        |u| + eps     2          |v| + eps     2  ]
+					 * J = --- \bar{v}^2 | log ( -----------  )   +  log ( -----------  )   |  
+					 *      2            [       |u    |+ eps              |v    |+ eps     ]
+					 *                              obs                       obs
+					 *        dJ                              1      u                             1
+					 * DU = - -- = - \bar{v}^2 log(u...) --------- ----  ~ - \bar{v}^2 log(u...) ------
+					 *        du                         |u| + eps  |u|                           u + eps
+					 */
+					for (i=0;i<NUMVERTICES;i++){
+						dux = - pow(meanvel,2.) * log((fabs(vx)+epsvel)/(fabs(vxobs)+epsvel)) / (vx+epsvel);
+						duy = - pow(meanvel,2.) * log((fabs(vy)+epsvel)/(fabs(vyobs)+epsvel)) / (vy+epsvel);
+						pe->values[i*NDOF4+0]+=dux*weight*Jdet*gauss->weight*l1l2l3[i]; 
+						pe->values[i*NDOF4+1]+=duy*weight*Jdet*gauss->weight*l1l2l3[i]; 
+					}
+					break;
+				default:
+					_error_("response %s not supported yet",EnumToStringx(responses[resp]));
+			}
 		}
 	}
@@ -4820,5 +4827,5 @@
 /*}}}*/
 /*FUNCTION Tria::SurfaceAbsVelMisfit {{{1*/
-double Tria::SurfaceAbsVelMisfit(bool process_units){
+double Tria::SurfaceAbsVelMisfit(bool process_units,int weight_index){
 
 	const int    numdof=NDOF2*NUMVERTICES;
@@ -4854,5 +4861,5 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -4910,5 +4917,5 @@
 /*}}}*/
 /*FUNCTION Tria::SurfaceAverageVelMisfit {{{1*/
-double Tria::SurfaceAverageVelMisfit(bool process_units){
+double Tria::SurfaceAverageVelMisfit(bool process_units,int weight_index){
 
 	const int    numdof=2*NUMVERTICES;
@@ -4945,5 +4952,5 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -4971,5 +4978,5 @@
 /*}}}*/
 /*FUNCTION Tria::SurfaceLogVelMisfit {{{1*/
-double Tria::SurfaceLogVelMisfit(bool process_units){
+double Tria::SurfaceLogVelMisfit(bool process_units,int weight_index){
 
 	const int    numdof=NDOF2*NUMVERTICES;
@@ -5008,5 +5015,5 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -5036,5 +5043,5 @@
 /*}}}*/
 /*FUNCTION Tria::SurfaceLogVxVyMisfit {{{1*/
-double Tria::SurfaceLogVxVyMisfit(bool process_units){
+double Tria::SurfaceLogVxVyMisfit(bool process_units,int weight_index){
 
 	const int    numdof=NDOF2*NUMVERTICES;
@@ -5073,5 +5080,5 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -5126,5 +5133,5 @@
 /*}}}*/
 /*FUNCTION Tria::SurfaceRelVelMisfit {{{1*/
-double Tria::SurfaceRelVelMisfit(bool process_units){
+double Tria::SurfaceRelVelMisfit(bool process_units,int weight_index){
 	const int  numdof=2*NUMVERTICES;
 
@@ -5162,5 +5169,5 @@
 
 		/*Get all parameters at gaussian point*/
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 		vx_input->GetParameterValue(&vx,gauss);
 		vy_input->GetParameterValue(&vy,gauss);
@@ -5190,5 +5197,5 @@
 /*}}}*/
 /*FUNCTION Tria::ThicknessAbsMisfit {{{1*/
-double Tria::ThicknessAbsMisfit(bool process_units){
+double Tria::ThicknessAbsMisfit(bool process_units,int weight_index){
 
 	/* Constants */
@@ -5230,5 +5237,5 @@
 		thickness_input->GetParameterDerivativeValue(&dH[0],&xyz_list[0][0],gauss);
 		thicknessobs_input->GetParameterValue(&thicknessobs,gauss);
-		weights_input->GetParameterValue(&weight,gauss,0);
+		weights_input->GetParameterValue(&weight,gauss,weight_index);
 
 		/*compute ThicknessAbsMisfit*/
Index: /issm/trunk/src/c/objects/Elements/Tria.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.h	(revision 8606)
+++ /issm/trunk/src/c/objects/Elements/Tria.h	(revision 8607)
@@ -123,10 +123,10 @@
 		void   MinVy(double* pminvy, bool process_units);
 		void   MinVz(double* pminvz, bool process_units);
-		double ThicknessAbsMisfit(bool process_units);
-		double SurfaceAbsVelMisfit(bool process_units);
-		double SurfaceRelVelMisfit(bool process_units);
-		double SurfaceLogVelMisfit(bool process_units);
-		double SurfaceLogVxVyMisfit(bool process_units);
-		double SurfaceAverageVelMisfit(bool process_units);
+		double ThicknessAbsMisfit(     bool process_units,int weight_index);
+		double SurfaceAbsVelMisfit(    bool process_units,int weight_index);
+		double SurfaceRelVelMisfit(    bool process_units,int weight_index);
+		double SurfaceLogVelMisfit(    bool process_units,int weight_index);
+		double SurfaceLogVxVyMisfit(   bool process_units,int weight_index);
+		double SurfaceAverageVelMisfit(bool process_units,int weight_index);
 		void   PatchFill(int* pcount, Patch* patch);
 		void   PatchSize(int* pnumrows, int* pnumvertices,int* pnumnodes);
Index: /issm/trunk/src/c/solutions/controltao_core.cpp
===================================================================
--- /issm/trunk/src/c/solutions/controltao_core.cpp	(revision 8606)
+++ /issm/trunk/src/c/solutions/controltao_core.cpp	(revision 8607)
@@ -113,5 +113,5 @@
 	//VecScale(gradient,-1.);
 	VecCopy(gradient,G);
-	CostFunctionx(fcn, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters,SurfaceAbsVelMisfitEnum);
+	CostFunctionx(fcn, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters);
 
 	//printf("X\n");
Index: /issm/trunk/src/c/solutions/objectivefunctionC.cpp
===================================================================
--- /issm/trunk/src/c/solutions/objectivefunctionC.cpp	(revision 8606)
+++ /issm/trunk/src/c/solutions/objectivefunctionC.cpp	(revision 8607)
@@ -25,12 +25,10 @@
 	
 	/*output: */
-	double J=0,Jplus;
+	double J;
 	
 	/*parameters: */
-	int        num_responses;
 	int        solution_type,analysis_type;
 	bool       isstokes       = false;
 	bool       conserve_loads = true;
-	int       *responses      = NULL;
 	FemModel  *femmodel       = NULL;
 
@@ -39,6 +37,4 @@
 
 	/*Recover parameters: */
-	femmodel->parameters->FindParam(&num_responses,NumResponsesEnum);
-	femmodel->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
 	femmodel->parameters->FindParam(&isstokes,IsStokesEnum);
 	femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum);
@@ -74,11 +70,7 @@
 
 	/*Compute misfit for this velocity field.*/
-	for(int i=0;i<num_responses;i++){
-		CostFunctionx(&Jplus, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters,responses[i]);
-		J+=Jplus;
-	}
+	CostFunctionx(&J, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters);
 
 	/*Free ressources:*/
-	xfree((void**)&responses);
 	return J;
 }
Index: /issm/trunk/src/m/solutions/objectivefunctionC.m
===================================================================
--- /issm/trunk/src/m/solutions/objectivefunctionC.m	(revision 8606)
+++ /issm/trunk/src/m/solutions/objectivefunctionC.m	(revision 8607)
@@ -6,6 +6,4 @@
 
 %recover some parameters
-num_responses = femmodel.parameters.NumResponses;
-responses     = femmodel.parameters.StepResponses;
 analysis_type = femmodel.parameters.AnalysisType;
 solution_type = femmodel.parameters.SolutionType;
@@ -36,5 +34,3 @@
 
 %Compute misfit for this velocity field
-for i=1:num_responses
-	J=J+CostFunction(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials, femmodel.parameters,responses(i));
-end
+J=CostFunction(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials, femmodel.parameters);
Index: /issm/trunk/src/mex/CostFunction/CostFunction.cpp
===================================================================
--- /issm/trunk/src/mex/CostFunction/CostFunction.cpp	(revision 8606)
+++ /issm/trunk/src/mex/CostFunction/CostFunction.cpp	(revision 8607)
@@ -14,5 +14,4 @@
 	Materials  *materials  = NULL;
 	Parameters *parameters = NULL;
-	int         response;
 
 	/* output datasets: */
@@ -32,5 +31,4 @@
 	FetchData((DataSet**)&materials,MATERIALS);
 	FetchParams(&parameters,PARAMETERS);
-	FetchData(&response,RESPONSE);
 
 	/*configure: */
@@ -40,5 +38,5 @@
 
 	/*!Call core code: */
-	CostFunctionx(&J, elements,nodes,vertices, loads,materials,parameters,response);
+	CostFunctionx(&J, elements,nodes,vertices, loads,materials,parameters);
 
 	/*write output : */
@@ -60,5 +58,5 @@
 {
 	_printf_(true,"\n");
-	_printf_(true,"   usage: [J] = %s(elements,nodes,vertices,loads, materials, parameters,response);\n",__FUNCT__);
+	_printf_(true,"   usage: [J] = %s(elements,nodes,vertices,loads, materials, parameters);\n",__FUNCT__);
 	_printf_(true,"\n");
 }
Index: /issm/trunk/src/mex/CostFunction/CostFunction.h
===================================================================
--- /issm/trunk/src/mex/CostFunction/CostFunction.h	(revision 8606)
+++ /issm/trunk/src/mex/CostFunction/CostFunction.h	(revision 8607)
@@ -23,5 +23,4 @@
 #define MATERIALS (mxArray*)prhs[4]
 #define PARAMETERS (mxArray*)prhs[5]
-#define RESPONSE (mxArray*)prhs[6]
 
 /* serial output macros: */
@@ -32,5 +31,5 @@
 #define NLHS  1
 #undef NRHS
-#define NRHS  7
+#define NRHS  6
 
 #endif  /* _COSTFUNCTION_H */
Index: /issm/trunk/src/mex/Response/Response.cpp
===================================================================
--- /issm/trunk/src/mex/Response/Response.cpp	(revision 8606)
+++ /issm/trunk/src/mex/Response/Response.cpp	(revision 8607)
@@ -16,4 +16,5 @@
 	char       *response   = NULL;
 	bool        process_units;
+	int         weight_index;
 
 	/* output datasets: */
@@ -35,4 +36,5 @@
 	FetchData(&response,RESPONSE);
 	FetchData(&process_units,PROCESSUNITS);
+	FetchData(&weight_index,WEIGHTINDEX);
 
 	/*configure: */
@@ -42,5 +44,5 @@
 
 	/*!Call core code: */
-	Responsex(&resp, elements,nodes,vertices, loads,materials,parameters,response,process_units);
+	Responsex(&resp, elements,nodes,vertices, loads,materials,parameters,response,process_units,weight_index);
 
 	/*write output : */
Index: /issm/trunk/src/mex/Response/Response.h
===================================================================
--- /issm/trunk/src/mex/Response/Response.h	(revision 8606)
+++ /issm/trunk/src/mex/Response/Response.h	(revision 8607)
@@ -26,4 +26,5 @@
 #define RESPONSE (mxArray*)prhs[6]
 #define PROCESSUNITS (mxArray*)prhs[7]
+#define WEIGHTINDEX (mxArray*)prhs[8]
 
 /* serial output macros: */
@@ -34,5 +35,5 @@
 #define NLHS  1
 #undef NRHS
-#define NRHS  8
+#define NRHS  9
 
 #endif  /* _RESPONSE_H */