From e924dca121b5b34bbbc4f1af76875ed1cc0ef97e Mon Sep 17 00:00:00 2001
From: Philippe Delandmeter <philippe.delandmeter@uclouvain.be>
Date: Fri, 10 May 2013 14:41:02 +0000
Subject: [PATCH] gmsh : new tranformation from a 2d field to LonLat
---
Mesh/Field.cpp | 53 ++++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 53 insertions(+)
diff --git a/Mesh/Field.cpp b/Mesh/Field.cpp
index b6ed27696d..39a400f01b 100644
--- a/Mesh/Field.cpp
+++ b/Mesh/Field.cpp
@@ -371,6 +371,58 @@ class LonLatField : public Field
}
};
+class XY2d2LonLatField : public Field
+{
+ int iField;
+ double R, phiP, thetaP;
+ public:
+ std::string getDescription()
+ {
+ return "Evaluate Field[IField] in geographic coordinates (longitude, latitude) from x and y projected on a plane\n\n";
+ }
+ XY2d2LonLatField()
+ {
+ iField = 1;
+ options["IField"] = new FieldOptionInt
+ (iField, "Index of the field to evaluate.");
+ R = 6371e3;
+ phiP = 0;
+ thetaP = 0;
+
+ options["Radius"] = new FieldOptionDouble
+ (R, "radius of the sphere");
+ options["Phi"] = new FieldOptionDouble
+ (phiP, "longitude of the projection point (in degrees)");
+ options["Theta"] = new FieldOptionDouble
+ (thetaP, "latitude of the projection point (in degrees)");
+ }
+ const char *getName()
+ {
+ return "XY2d2LonLat";
+ }
+ double operator() (double x2d, double y2d, double z2d, GEntity *ge=0)
+ {
+ Field *field = GModel::current()->getFields()->get(iField);
+ if(!field || iField == id) return MAX_LC;
+ double phi = phiP*M_PI/180;
+ double theta = thetaP*M_PI/180;
+ double pOx = cos(theta)*cos(phi)*R;
+ double pOy = cos(theta)*sin(phi)*R;
+ double pOz = sin(theta)*R;
+ double pPhiX = -sin(phi);
+ double pPhiY = cos(phi);
+ double pPhiZ = 0;
+ double pThetaX = -sin(theta)*cos(phi);
+ double pThetaY = -sin(theta)*sin(phi);
+ double pThetaZ = cos(theta);
+
+ double x = pPhiX * x2d + pThetaX * y2d + pOx;
+ double y = pPhiY * x2d + pThetaY * y2d + pOy;
+ double z = pPhiZ * x2d + pThetaZ * y2d + pOz;
+ return (*field)(atan2(y, x), asin(z / R), 0);
+ }
+};
+
class BoxField : public Field
{
double v_in, v_out, x_min, x_max, y_min, y_max, z_min, z_max;
@@ -2063,6 +2115,7 @@ FieldManager::FieldManager()
map_type_name["Cylinder"] = new FieldFactoryT<CylinderField>();
map_type_name["Frustum"] = new FieldFactoryT<FrustumField>();
map_type_name["LonLat"] = new FieldFactoryT<LonLatField>();
+ map_type_name["XY2d2LonLat"] = new FieldFactoryT<XY2d2LonLatField>();
#if defined(HAVE_POST)
map_type_name["PostView"] = new FieldFactoryT<PostViewField>();
#endif
--
GitLab