Changeset 12600
- Timestamp:
- 07/03/12 08:39:17 (13 years ago)
- Location:
- issm/trunk-jpl/src/c
- Files:
-
- 9 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.cpp
r12511 r12600 12 12 using namespace std; 13 13 14 int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,double* index,double* x,double* y,int nods,int nels){ 15 16 /*Intermediary*/ 17 int i,j,k; 18 int verbose=0; 19 int noerr=1; 14 int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,int* index,double* x,double* y,int nods,int nels){ 20 15 21 16 /*Options*/ 22 BamgOpts* bamgopts=NULL; 23 bamgopts=new BamgOpts(); 17 BamgOpts* bamgopts=new BamgOpts(); 24 18 25 // read mesh 26 if(verbose) _printLine_("Reading mesh"); 19 /*read mesh*/ 27 20 Mesh Th(index,x,y,nods,nels); 28 21 29 //write mesh and geometry 30 if (verbose) _printLine_("Write Geometry"); 22 /*write mesh and geometry*/ 31 23 Th.Gh.WriteGeometry(bamggeom,bamgopts); 32 if (verbose) _printLine_("Write Mesh");33 24 Th.WriteMesh(bamgmesh,bamgopts); 34 25 35 / /clean up26 /*clean up and return*/ 36 27 delete bamgopts; 37 38 /*No error return*/ 39 return noerr; 28 return 1; 40 29 41 30 } -
issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.h
r3913 r12600 9 9 10 10 /* local prototypes: */ 11 int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom, double* index,double* x,double* y,int nods,int nels);11 int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,int* index,double* x,double* y,int nods,int nels); 12 12 13 13 #endif -
issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp
r12511 r12600 13 13 using namespace std; 14 14 15 int InterpFromMeshToMesh2dx(double** pdata_interp, double* index_data,double* x_data,double* y_data,int nods_data,int nels_data,16 double* data,int data_rows,int data_cols,double* x_interp,double* y_interp,int nods_interp,double* default_values,int num_default_values, DataSet* contours){15 int InterpFromMeshToMesh2dx(double** pdata_interp,int* index_data,double* x_data,double* y_data,int nods_data,int nels_data, 16 double* data,int M_data,int N_data,double* x_interp,double* y_interp,int N_interp,Options* options){ 17 17 18 18 /*Output*/ … … 20 20 21 21 /*Intermediary*/ 22 bool isdefault; 23 double defaultvalue; 22 24 R2 r; 23 25 I2 I; … … 29 31 double data_value; 30 32 Icoor2 dete[3]; 31 int verbose=0;32 33 /*default values: */34 Vector* vec_incontour=NULL;35 double* incontour=NULL;36 bool skip_bamg=false;37 33 38 34 /*Checks*/ 39 if ( data_cols<=0){40 _error2_("data provided has a negative number of columns");35 if (M_data!=nods_data && M_data!=nels_data){ 36 _error2_("data provided should have either " << nods_data << " or " << nels_data << " lines (not " << M_data << ")"); 41 37 } 42 if (data_rows!=nods_data && data_rows!=nels_data){ 43 _error2_("data provided should have either " << nods_data << " or " << nels_data << " lines (not " << data_rows << ")"); 38 39 /*Get default*/ 40 if(options->GetOption("default")){ 41 isdefault=true; 42 options->Get(&defaultvalue,"default"); 44 43 } 45 if((num_default_values) && (data_cols>1)){ 46 _error2_("data provided can only have 1 column if a default value is provided"); 47 } 48 49 /*If default values supplied, figure out which nodes are inside the contour, including the border of the contour: */ 50 if(num_default_values){ 51 ContourToNodesx( &vec_incontour,x_interp,y_interp,nods_interp,contours,1); 52 incontour=vec_incontour->ToMPISerial(); 44 else{ 45 isdefault=false; 53 46 } 54 47 55 48 /*Initialize output*/ 56 if (verbose) _printLine_("Initializing output vector"); 57 data_interp=xNew<double>(nods_interp*data_cols); 49 data_interp=xNew<double>(N_interp*N_data); 58 50 59 // read background mesh 60 if (verbose) _printLine_("Reading mesh"); 51 /*read background mesh*/ 61 52 Mesh Th(index_data,x_data,y_data,nods_data,nels_data); 62 53 Th.CreateSingleVertexToTriangleConnectivity(); 63 54 64 //Loop over output nodes 65 if (verbose) _printLine_("Loop over the nodes"); 66 for(i=0;i<nods_interp;i++){ 67 68 /*reset skip_bamg: */ 69 skip_bamg=false; 55 /*Loop over output nodes*/ 56 for(i=0;i<N_interp;i++){ 70 57 71 /*figure out if we should skip bamg logic: */ 72 if(num_default_values){ 73 if(!incontour[i]){ 74 /*This node is not inside the contour. Skip Bamg logic and apply default value.: */ 75 skip_bamg=true; 58 //Get current point coordinates 59 r.x=x_interp[i]; r.y=y_interp[i]; 60 I2 I=Th.R2ToI2(r); 61 62 //Find triangle holding r/I 63 Triangle &tb=*Th.TriangleFindFromCoord(I,dete); 64 65 // internal point 66 if (tb.det>0){ 67 //Area coordinate 68 areacoord[0]= (double) dete[0]/tb.det; 69 areacoord[1]= (double) dete[1]/tb.det; 70 areacoord[2]= (double) dete[2]/tb.det; 71 //3 vertices of the triangle 72 i0=Th.GetId(tb[0]); 73 i1=Th.GetId(tb[1]); 74 i2=Th.GetId(tb[2]); 75 //triangle number 76 it=Th.GetId(tb); 77 } 78 //external point 79 else{ 80 if(isdefault){ 81 for(j=0;j<N_data;j++) data_interp[i*N_data+j]=defaultvalue; 82 continue; 76 83 } 77 } 78 79 if(skip_bamg==false){ 80 81 //Get current point coordinates 82 r.x=x_interp[i]; r.y=y_interp[i]; 83 I2 I=Th.R2ToI2(r); 84 85 //Find triangle holding r/I 86 Triangle &tb=*Th.TriangleFindFromCoord(I,dete); 87 88 // internal point 89 if (tb.det>0){ 90 //Area coordinate 91 areacoord[0]= (double) dete[0]/tb.det; 92 areacoord[1]= (double) dete[1]/tb.det; 93 areacoord[2]= (double) dete[2]/tb.det; 94 //3 vertices of the triangle 95 i0=Th.GetId(tb[0]); 96 i1=Th.GetId(tb[1]); 97 i2=Th.GetId(tb[2]); 98 //triangle number 99 it=Th.GetId(tb); 100 } 101 //external point 102 else { 84 else{ 103 85 //Get closest adjacent triangle (inside the mesh) 104 86 AdjacentTriangle ta=CloseBoundaryEdge(I,&tb,aa,bb).Adj(); … … 116 98 it=Th.GetId(tc); 117 99 } 118 119 if (data_rows==nods_data){ 120 for (j=0;j<data_cols;j++){ 121 data_interp[i*data_cols+j]=areacoord[0]*data[data_cols*i0+j]+areacoord[1]*data[data_cols*i1+j]+areacoord[2]*data[data_cols*i2+j]; 122 } 123 } 124 else{ 125 for (j=0;j<data_cols;j++){ 126 if (it<0 || it>=nels_data){ 127 _error2_("Triangle number " << it << " not in [0 " << nels_data << "], because not correctly implemented yet... interpolate on grid first"); 128 } 129 data_interp[i*data_cols+j]=data[data_cols*it+j]; 130 } 100 } 101 102 if (M_data==nods_data){ 103 for (j=0;j<N_data;j++){ 104 data_interp[i*N_data+j]=areacoord[0]*data[N_data*i0+j]+areacoord[1]*data[N_data*i1+j]+areacoord[2]*data[N_data*i2+j]; 131 105 } 132 106 } 133 107 else{ 134 if(num_default_values==1) data_interp[i]=default_values[0]; 135 else data_interp[i]=default_values[i]; 108 for (j=0;j<N_data;j++){ 109 if (it<0 || it>=nels_data){ 110 _error2_("Triangle number " << it << " not in [0 " << nels_data << "], report bug to developers"); 111 } 112 data_interp[i*N_data+j]=data[N_data*it+j]; 113 } 136 114 } 137 115 } 138 116 139 117 /*Assign output pointers:*/ 140 if (verbose) _printLine_("Assigning output");141 118 *pdata_interp=data_interp; 142 119 -
issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.h
r12121 r12600 6 6 #define _INTERPFROMMESHTOMESH2DX_H 7 7 8 #include "../../objects/objects.h" 8 class Options; 9 9 10 /* local prototypes: */ 11 int InterpFromMeshToMesh2dx(double** pdata_interp,double* index_data,double* x_data,double* y_data,int nods_data,int nels_data, 12 double* data,int data_rows,int data_cols,double* x_interp,double* y_interp,int nods_interp,double* default_values,int num_default_values,DataSet* contours); 10 int InterpFromMeshToMesh2dx(double** pdata_interp,int* index_data,double* x_data,double* y_data,int nods_data,int nels_data, 11 double* data,int M_data,int N_data,double* x_interp,double* y_interp,int N_interp,Options* options); 13 12 14 13 #endif -
issm/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp
r12595 r12600 62 62 double *data = NULL; 63 63 int *index = NULL; 64 double *indexd= NULL;65 64 66 65 observations->ObservationList(&x,&y,&data,&nobs); … … 68 67 _printLine_("Generation Delaunay Triangulation"); 69 68 BamgTriangulatex(&index,&nel,x,y,nobs); 70 indexd =xNewZeroInit<double>(nel*3);71 for(int i=0;i<nel*3;i++) indexd[i]=(double)index[i];72 xDelete<int>(index);73 69 74 70 _printLine_("Interpolating"); 75 71 xDelete<double>(predictions); 76 InterpFromMeshToMesh2dx(&predictions,index d,x,y,nobs,nel,data,nobs,1,x_interp,y_interp,n_interp,NULL,0,new DataSet());72 InterpFromMeshToMesh2dx(&predictions,index,x,y,nobs,nel,data,nobs,1,x_interp,y_interp,n_interp,options); 77 73 xDelete<double>(x); 78 74 xDelete<double>(y); -
issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.cpp
r12450 r12600 13 13 using namespace std; 14 14 15 void TriaSearchx(double** ptria, double* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes){15 void TriaSearchx(double** ptria,int* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes){ 16 16 17 17 /*Output*/ -
issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.h
r8303 r12600 9 9 10 10 /* local prototypes: */ 11 void TriaSearchx(double** ptria, double* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes);11 void TriaSearchx(double** ptria,int* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes); 12 12 13 13 #endif -
issm/trunk-jpl/src/c/objects/Bamg/Mesh.cpp
r12574 r12600 41 41 } 42 42 /*}}}*/ 43 /*FUNCTION Mesh::Mesh( double* index,double* x,double* y,int nods,int nels){{{*/44 Mesh::Mesh( double* index,double* x,double* y,int nods,int nels):Gh(*(new Geometry())),BTh(*this){43 /*FUNCTION Mesh::Mesh(int* index,double* x,double* y,int nods,int nels){{{*/ 44 Mesh::Mesh(int* index,double* x,double* y,int nods,int nels):Gh(*(new Geometry())),BTh(*this){ 45 45 46 46 Init(0); … … 265 265 266 266 /*IO*/ 267 /*FUNCTION Mesh::ReadMesh( double* index,double* x,double* y,int nods,int nels){{{*/268 void Mesh::ReadMesh( double* index,double* x,double* y,int nods,int nels){267 /*FUNCTION Mesh::ReadMesh(int* index,double* x,double* y,int nods,int nels){{{*/ 268 void Mesh::ReadMesh(int* index,double* x,double* y,int nods,int nels){ 269 269 270 270 double Hmin = HUGE_VAL;// the infinie value -
issm/trunk-jpl/src/c/objects/Bamg/Mesh.h
r12218 r12600 57 57 //Constructors/Destructors 58 58 Mesh(BamgGeom* bamggeom,BamgMesh* bamgmesh,BamgOpts* bamgopts); 59 Mesh( double* index,double* x,double* y,int nods,int nels);/*MeshConvert*/59 Mesh(int* index,double* x,double* y,int nods,int nels);/*MeshConvert*/ 60 60 Mesh(double* x,double* y,int nods); /*BamgTriangulate*/ 61 61 Mesh(Mesh &,Geometry * pGh=0,Mesh* pBTh=0,long maxnbv_in=0 ); //copy operator … … 110 110 BamgVertex* NearestVertex(Icoor1 i,Icoor1 j) ; 111 111 Triangle* TriangleFindFromCoord(const I2 & ,Icoor2 [3],Triangle *tstart=0) const; 112 void ReadMesh( double* index,double* x,double* y,int nods,int nels);112 void ReadMesh(int* index,double* x,double* y,int nods,int nels); 113 113 void ReadMesh(BamgMesh* bamgmesh, BamgOpts* bamgopts); 114 114 void WriteMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts);
Note:
See TracChangeset
for help on using the changeset viewer.