Salome HOME
Correction for hydro_test
[modules/hydro.git] / src / HYDROTools / shapesGroups.py
index 18dedea56ae3a3a7b25330e9142302fa090a6a44..7f3b7e0b99794aaa46ea6c239253c48f092a2746 100644 (file)
@@ -12,7 +12,7 @@ import medcoupling as mc
 import shapefile
 import math
 import os
-
+import json
 
 def freeBordersGroup(meshFileIn, meshFileOut=""):
     """
@@ -23,7 +23,8 @@ def freeBordersGroup(meshFileIn, meshFileOut=""):
     return full path of the output mesh file
     """
     smesh = smeshBuilder.New()
-    smesh.SetEnablePublish( False ) # Set to False to avoid publish in study if not needed
+    if not salome.sg.hasDesktop():
+        smesh.SetEnablePublish( False ) # Set to False to avoid publish in study if not needed
     ([MESH], status) = smesh.CreateMeshesFromMED(meshFileIn)
 
     nbAdded, MESH, addedBnd = MESH.MakeBoundaryElements( SMESH.BND_1DFROM2D, '', '', 0, [])
@@ -158,6 +159,7 @@ def writeShapeLines(mcMesh, grpName, nodeChains, outputDirectory, offsetX=0., of
         w.line([chaincoords])
         w.record(chainName)
     w.close()
+    return shapeFileName + '.shp'
 
 
 def writeShapePoints(mcMesh, grpName, nodeChains, outputDirectory, offsetX=0., offsetY=0.):
@@ -206,7 +208,7 @@ def exploreEdgeGroups(meshFile, outputDirectory="", offsetX=0., offsetY=0.):
     prefix = os.path.basename(a[0]) # prefix = file name without extension
 
     mcMesh = ml.MEDFileMesh.New(meshFile)
-    dim = mcMesh.getSpaceDimension()
+    dim = mcMesh.getMeshDimension()
     d1=-1        # when dimension 2, edges are dim -1
     if dim == 3: # when dimension 3, edges are dim -2
         d1=-2
@@ -215,20 +217,25 @@ def exploreEdgeGroups(meshFile, outputDirectory="", offsetX=0., offsetY=0.):
 
     groups = [mcMesh.getGroup(d1, name) for name in grp_names] # list of groups in their name order
 
+    filenames = []
     for (grp, grpName) in zip(groups, grp_names):
         fullGrpName = prefix + '_' + grpName
         nodeChains = explodeGroup(grp, fullGrpName)
-        writeShapeLines(mcMesh, fullGrpName, nodeChains, outputDirectory, offsetX, offsetY)
+        filename = writeShapeLines(mcMesh, fullGrpName, nodeChains, outputDirectory, offsetX, offsetY)
         writeShapePoints(mcMesh, fullGrpName, nodeChains, outputDirectory, offsetX, offsetY)
+        filenames.append(filename)
+    shapesListFile = os.path.join(outputDirectory, "shapesList.json")
+    with open(shapesListFile, 'w') as f:
+        json.dump(filenames, f)
 
-
-def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder=True, splitShapeAdjusted=True):
+def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, outputDirectory="", splitFreeBorder=True, splitShapeAdjusted=True):
     """
     shapeFileToAdjust must be a closed line or polygon crossing freeBorderShapefile in 2 points.
     Find in shapeFileToAdjust and in freeBorderShapefile the two closest corresponding points and move the points in shapeFileToAdjust to correspond to the points found in freeBorderShapefile. Split shapeFileToAdjust in two parts (inside or outside freeBorder). If requested, split freeBorderShapefile in two parts. Same for shapeFileToAdjust.
     parameters:
     freeBorderShapefile: a set of free border lines, as generated by the functions freeBordersGroup and exploreEdgeGroups.
     shapefileToAdjust: a closed line or polygon, supposed to be drawn in qgis to pass as close as possible to the points to be connected, on the free border.
+    outputDirectory: if empty, write the resulting shapefiles in their respective directory, with the suffix '_adj', otherwise write in the outputDirectory.
     splitFreeBorder: boolean default True
     splitShapeAdjusted: boolean default True
     """
@@ -257,12 +264,23 @@ def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder
     pdist =[]
     x0 = -1.e30
     y0 = -1.e30
+    discountLastSftaPoint = 0
     for i,p in enumerate(sfta.points):
-        d = math.sqrt((x0-p[0])*(x0-p[0]) + (y0-p[1])*(y0-p[1]))
-        x0 = p[0]
+        if i == 0:
+            x00 = p[0] # keep first point
+            y00 = p[1]
+        d = math.sqrt((x0-p[0])*(x0-p[0]) + (y0-p[1])*(y0-p[1])) # distance to previous point
+        x0 = p[0] # keep previous point
         y0 = p[1]
         if d < 1.e-5:
-            continue # do not take into account superposed points in shape to adjust
+            print("two consecutives points of shapefile To adjust are superposed, OK")
+            continue # do not take into account consecutive superposed points in shape to adjust
+        if i == len(sfta.points) -1:
+            d = math.sqrt((x00-p[0])*(x00-p[0]) + (y00-p[1])*(y00-p[1])) # distance between first and last point
+            if d < 1.e-5:
+                discountLastSftaPoint = 1
+                print("last point of shapefile To adjust is superposed to first point, last point discarded, OK")
+                continue # do not take into account last point if superposed to first point, in shape to adjust 
         dmin = 1.e30
         jmin = -1
         for j,pfb in enumerate(fbs.points):
@@ -287,9 +305,11 @@ def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder
     # --- write the adusted shapefile: free border closest points replaced with corresponding points
     #     on the free border. two polylines, one inside the domain, one outside
 
-    a = os.path.splitext(shapefileToAdjust)
-    shapefileAdjusted = a[0] + '_adj' + a[1]
-    chainName = os.path.basename(a[0])
+    if outputDirectory == "":
+        outputDirectory = os.path.dirname(shapefileToAdjust)
+    a = os.path.splitext(os.path.basename(shapefileToAdjust))
+    shapefileAdjusted = os.path.join(outputDirectory, a[0] + '_adj' + a[1])
+    chainName = a[0] + '_adj'
 
     w = shapefile.Writer(shapefileAdjusted)
     w.shapeType = 3
@@ -307,7 +327,7 @@ def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder
         chaincoords = []
         chaincoords.append(fbs.points[ifb2])
         if i2+1 < len(sfta.points):
-            for i in range(i2+1, len(sfta.points)):
+            for i in range(i2+1, len(sfta.points) -discountLastSftaPoint):
                 chaincoords.append(sfta.points[i])
         for i in range(i1):
             chaincoords.append(sfta.points[i])
@@ -321,10 +341,12 @@ def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder
             chaincoords.append(sfta.points[i])
         chaincoords.append(fbs.points[ifb2])
         if i2+1 < len(sfta.points):
-            for i in range(i2+1, len(sfta.points)):
+            for i in range(i2+1, len(sfta.points) -discountLastSftaPoint):
                 chaincoords.append(sfta.points[i])
         for i in range(i1):
             chaincoords.append(sfta.points[i])
+        if discountLastSftaPoint:
+            chaincoords.append(fbs.points[ifb1]) # close shape when first point if superposed with last
         w.line([chaincoords])
         w.record(chainName)
 
@@ -333,9 +355,11 @@ def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust, splitFreeBorder
     if splitFreeBorder:
         # write the free border split in two polylines (cut by the adjusted shapefile)
 
-        a = os.path.splitext(freeBorderShapefile)
-        freeBorderSplit = a[0] + '_split' + a[1]
-        chainName = os.path.basename(a[0])
+        if outputDirectory == "":
+            outputDirectory = os.path.dirname(freeBorderShapefile)
+        a = os.path.splitext(os.path.basename(freeBorderShapefile))
+        freeBorderSplit = os.path.join(outputDirectory, a[0] + '_split' + a[1])
+        chainName = a[0] + '_split'
 
         w = shapefile.Writer(freeBorderSplit)
         w.shapeType = 3