//  SMESH DriverMED : driver to read and write 'med' files
//
//  Copyright (C) 2003  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
//  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS 
// 
//  This library is free software; you can redistribute it and/or 
//  modify it under the terms of the GNU Lesser General Public 
//  License as published by the Free Software Foundation; either 
//  version 2.1 of the License. 
// 
//  This library is distributed in the hope that it will be useful, 
//  but WITHOUT ANY WARRANTY; without even the implied warranty of 
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU 
//  Lesser General Public License for more details. 
// 
//  You should have received a copy of the GNU Lesser General Public 
//  License along with this library; if not, write to the Free Software 
//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA 
// 
//  See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org 
//
//
//
//  File   : DriverMED_W_SMESHDS_Mesh.cxx
//  Module : SMESH

using namespace std;
#include "DriverMED_W_SMESHDS_Mesh.h"
#include "DriverMED_W_SMDS_Mesh.h"

#include "SMDS_MeshElement.hxx"
#include "SMDS_MeshNode.hxx"
#include "SMDS_MeshFacesIterator.hxx"
#include "SMDS_MeshEdgesIterator.hxx"
#include "SMDS_MeshNodesIterator.hxx"
#include "SMDS_MeshVolumesIterator.hxx"
#include "SMESHDS_SubMesh.hxx"
#include "SMDS_MapIteratorOfExtendedMap.hxx"
#include "SMDS_MapOfMeshElement.hxx"
#include "SMDS_Position.hxx"

#include <TopExp.hxx>
#include <TopExp_Explorer.hxx>
#include <TopoDS.hxx>
#include <TopoDS_Iterator.hxx>
#include <TopoDS_Compound.hxx>
#include <TopoDS_CompSolid.hxx>
#include <TopoDS_Solid.hxx>
#include <TopoDS_Shell.hxx>
#include <TopoDS_Face.hxx>
#include <TopoDS_Wire.hxx>
#include <TopoDS_Edge.hxx>
#include <TopoDS_Vertex.hxx>
#include <TopoDS_Shape.hxx>
#include <TopTools_MapOfShape.hxx>
#include <TColStd_ListIteratorOfListOfInteger.hxx>


#include <map>
#include <set>
#include <vector>
#include "utilities.h"

DriverMED_W_SMESHDS_Mesh::DriverMED_W_SMESHDS_Mesh() {
  ;
}

DriverMED_W_SMESHDS_Mesh::~DriverMED_W_SMESHDS_Mesh() {
  ;
}

void DriverMED_W_SMESHDS_Mesh::SetMesh(Handle(SMDS_Mesh)& aMesh) {
  myMesh = aMesh;
}

void DriverMED_W_SMESHDS_Mesh::SetFile(string aFile) {
  myFile = aFile;
}

void DriverMED_W_SMESHDS_Mesh::SetFileId(med_idt aFileId) {
  myFileId = aFileId;
}

void DriverMED_W_SMESHDS_Mesh::SetMeshId(int aMeshId) {
  myMeshId = aMeshId;
}

void DriverMED_W_SMESHDS_Mesh::Write() {

  string myClass = string("SMDS_Mesh");
  string myExtension = string("MED");

  DriverMED_W_SMDS_Mesh* myWriter = new DriverMED_W_SMDS_Mesh;

  myWriter->SetMesh(myMesh);
  //  myWriter->SetFile(myFile);
  myWriter->SetMeshId(myMeshId);
  myWriter->SetFileId(myFileId);

  myWriter->Write();


}

void DriverMED_W_SMESHDS_Mesh::Add() {

  med_err ret = 0;
  int i,j,k,l;
  int numero;
  char message[200];
  Standard_Boolean ok;
  /* nombre d'objets MED */
  char nom_universel[MED_TAILLE_LNOM+1];
  med_int long_fichier_en_tete;
  char *fichier_en_tete;
  char version_hdf[10];
  char version_med[10];
  med_int nmaa,mdim,nnoe;
  med_int nmai[MED_NBR_GEOMETRIE_MAILLE],nfac[MED_NBR_GEOMETRIE_FACE];
  med_int nare[MED_NBR_GEOMETRIE_ARETE];
  /* nom du maillage */
  char nommaa[MED_TAILLE_NOM+1];
  /* noeuds */
  med_float *coo;
  // PN : Initilialisation de nomcoo et unicoo pour lisibilite du maillage
  char nomcoo[3*MED_TAILLE_PNOM+1]="x       y        z      ";
  char unicoo[3*MED_TAILLE_PNOM+1]="m       m        m     ";
  char *nomnoe;
  med_int *numnoe;
  med_int *nufano;
  med_repere rep;
  med_booleen inonoe,inunoe;
  med_mode_switch mode_coo;
  char str[MED_TAILLE_PNOM+1];
  med_int nbNodes;
  /* elements */
  med_int nsup;
  med_int edim;
  med_int taille;
  med_int elem_id,myId;
  med_int *connectivite;
  char *nomele;
  med_int *numele;
  med_int *nufael;
  med_booleen inoele, inuele;
  med_connectivite typ_con;
  med_geometrie_element typgeo;
  med_geometrie_element typmai[MED_NBR_GEOMETRIE_MAILLE] = {MED_POINT1,MED_SEG2,
							    MED_SEG3,MED_TRIA3,
							    MED_TRIA6,MED_QUAD4,
							    MED_QUAD8,MED_TETRA4,
							    MED_TETRA10,MED_HEXA8,
							    MED_HEXA20,MED_PENTA6,
							    MED_PENTA15,MED_PYRA5,
							    MED_PYRA13};
  med_int desmai[MED_NBR_GEOMETRIE_MAILLE] = {0,2,3,3,3,4,4,4,4,6,6,5,5,5,5};
  med_int nmailles[MED_NBR_GEOMETRIE_MAILLE];
  char nommai[MED_NBR_GEOMETRIE_MAILLE] [MED_TAILLE_NOM+1] = {"MED_POINT1",
							      "MED_SEG2",
							      "MED_SEG3",
							      "MED_TRIA3",
							      "MED_TRIA6",
							      "MED_QUAD4",
							      "MED_QUAD8",
							      "MED_TETRA4",
							      "MED_TETRA10",
							      "MED_HEXA8",
							      "MED_HEXA20",
							      "MED_PENTA6",
							      "MED_PENTA15",
							      "MED_PYRA5",
							      "MED_PYRA13"};
  med_geometrie_element typfac[MED_NBR_GEOMETRIE_FACE] = {MED_TRIA3,MED_TRIA6,
							  MED_QUAD4,MED_QUAD8};
  med_int desfac[MED_NBR_GEOMETRIE_FACE] = {3,3,4,4};
  med_int nfaces[MED_NBR_GEOMETRIE_FACE];
  char nomfac[MED_NBR_GEOMETRIE_FACE][MED_TAILLE_NOM+1] = {"MED_TRIA3","MED_TRIA6",
							   "MED_QUAD4","MED_QUAD8"};
  med_geometrie_element typare[MED_NBR_GEOMETRIE_ARETE] = {MED_SEG2,MED_SEG3};
  med_int desare[MED_NBR_GEOMETRIE_ARETE] = {2,3};
  med_int naretes[MED_NBR_GEOMETRIE_ARETE];
  char nomare[MED_NBR_GEOMETRIE_ARETE] [MED_TAILLE_NOM+1] = {"MED_SEG2","MED_SEG3"};

  typ_con = MED_NOD;
  mode_coo = MED_FULL_INTERLACE;
  numero = myMeshId;

  //---- provisoire : switch pour ecrire les familles de mailles
  int besoinfamilledemaille=1;
  //---- provisoire : switch pour ecrire les familles de mailles

  /****************************************************************************
   *                      OUVERTURE DU FICHIER EN ECRITURE                    *
   ****************************************************************************/
  char* file2Read = (char*)myFile.c_str();

  MESSAGE ( " file2Read " << file2Read )

    myFileId = MEDouvrir(file2Read,MED_REMP);
  if (myFileId < 0)
    {
      fprintf(stderr,">> ERREUR : ouverture du fichier %s \n",file2Read);
      exit(EXIT_FAILURE);
    }

  /****************************************************************************
   *                       NOMBRES D'OBJETS MED                               *
   ****************************************************************************/
  MESSAGE("(****************************)");
  MESSAGE("(* INFORMATIONS GENERALES : *)");
  MESSAGE("(****************************)");

  /* calcul de la dimension */
  mdim=2;
  double epsilon=0.00001;
  double nodeRefX;
  double nodeRefY;
  double nodeRefZ;

  bool dimX = true;
  bool dimY = true;
  bool dimZ = true;

  SMDS_MeshNodesIterator myItNodes(myMesh);
  int inode = 0;
  for (;myItNodes.More();myItNodes.Next()) {
    const Handle(SMDS_MeshElement)& elem = myItNodes.Value();
    const Handle(SMDS_MeshNode)& node = myMesh->GetNode(1,elem);
    if ( inode == 0 ) {
      nodeRefX = fabs(node->X());
      nodeRefY = fabs(node->Y());
      nodeRefZ = fabs(node->Z());
    }
    SCRUTE( inode );
    SCRUTE( nodeRefX );
    SCRUTE( nodeRefY );
    SCRUTE( nodeRefZ );

    if ( inode !=0 ) {
      if ( (fabs(fabs(node->X()) - nodeRefX) > epsilon ) && dimX )
	dimX = false;
      if ( (fabs(fabs(node->Y()) - nodeRefY) > epsilon ) && dimY )
	dimY = false;
      if ( (fabs(fabs(node->Z()) - nodeRefZ) > epsilon ) && dimZ )
	dimZ = false;
    }
    if ( !dimX && !dimY && !dimZ ) {
      mdim = 3;
      break;
    }
    inode++;
  }

  if ( mdim != 3 ) {
    if ( dimX && dimY && dimZ )
      mdim = 0;
    else if ( !dimX ) {
      if ( dimY && dimZ )
	mdim = 1;
      else if (( dimY && !dimZ ) || ( !dimY && dimZ ) )
	mdim = 2;
    } else if ( !dimY ) {
      if ( dimX && dimZ )
	mdim = 1;
      else if (( dimX && !dimZ ) || ( !dimX && dimZ ) )
	mdim = 2;
    } else if ( !dimZ ) {
      if ( dimY && dimX )
	mdim = 1;
      else if (( dimY && !dimX ) || ( !dimY && dimX ) )
	mdim = 2;
    }
  }

  MESSAGE ( " mdim " << mdim );

  /* creation du maillage */
  //mdim=3;
  sprintf(nommaa,"Mesh %d",numero);
  SCRUTE(nommaa);
  ret = MEDmaaCr(myFileId,nommaa,mdim);

  ASSERT(ret == 0);
  SCRUTE( ret );

  /* Combien de noeuds ? */
  nnoe = myMesh->NbNodes();
  //SCRUTE(nnoe);
  /* Combien de mailles, faces ou aretes ? */
  for (i=0;i<MED_NBR_GEOMETRIE_MAILLE;i++)
    nmailles[i]=0;

  Standard_Integer nb_of_nodes, nb_of_faces, nb_of_edges;
  vector<int> elem_Id[MED_NBR_GEOMETRIE_MAILLE];

  SMDS_MeshEdgesIterator itEdges(myMesh);
  nb_of_edges = myMesh->NbEdges();
  for (;itEdges.More();itEdges.Next()) {
    const Handle(SMDS_MeshElement)& elem = itEdges.Value();

    nb_of_nodes = elem->NbNodes();

    switch (nb_of_nodes) {
    case 2 : {
      elem_Id[1].push_back(elem->GetID());
      nmailles[1]++;
      break;
    }
    case 3 : {
      elem_Id[2].push_back(elem->GetID());
      nmailles[2]++;
      break;
    }
    }
  }


  SMDS_MeshFacesIterator itFaces(myMesh);
  nb_of_faces = myMesh->NbFaces();
  for (;itFaces.More();itFaces.Next()) {
    const Handle(SMDS_MeshElement)& elem = itFaces.Value();

    nb_of_nodes = elem->NbNodes();

    switch (nb_of_nodes) {
    case 3 : {
      elem_Id[3].push_back(elem->GetID());
      nmailles[3]++;
      break;
    }
    case 4 : {
      elem_Id[5].push_back(elem->GetID());
      nmailles[5]++;
      break;
    }
    case 6 : {
      elem_Id[4].push_back(elem->GetID());
      nmailles[4]++;
      break;
    }
    }

  }

  SMDS_MeshVolumesIterator itVolumes(myMesh);
  for (;itVolumes.More();itVolumes.Next()) {
    const Handle(SMDS_MeshElement)& elem = itVolumes.Value();

    nb_of_nodes = elem->NbNodes();
    switch (nb_of_nodes) {
    case 8 : {
      elem_Id[9].push_back(elem->GetID());
      nmailles[9]++;
      break;
    }
    }
  }



  /****************************************************************************
   *                       ECRITURE DES NOEUDS                                *
   ****************************************************************************/
  MESSAGE("(************************)");
  MESSAGE("(* NOEUDS DU MAILLAGE : *)");
  MESSAGE("(************************)");

  /* Allocations memoires */
  /* table des coordonnees
     profil : (dimension * nombre de noeuds ) */
  coo = (med_float*) malloc(sizeof(med_float)*nnoe*mdim);
  /* table  des numeros, des numeros de familles des noeuds
     profil : (nombre de noeuds) */
  numnoe = (med_int*) malloc(sizeof(med_int)*nnoe);
  nufano = (med_int*) malloc(sizeof(med_int)*nnoe);
  /* table des noms des noeuds
     profil : (nnoe*MED_TAILLE_PNOM+1) */
  nomnoe ="";

  /* PN  pour aster, il faut une famille 0 pour les noeuds et une autre pour les elements*/
  /* PN : Creation de la famille 0 */
  char * nomfam="FAMILLE_0";
  char *attdes="";
  char *gro=0;
  med_int ngro=0;
  med_int natt=1;
  med_int attide = 0;
  med_int attval=0;
  med_int numfam = 0;
  med_int attvalabs=1;
  ret = MEDfamCr(myFileId,nommaa,nomfam,numfam, &attide,&attval,attdes,natt,gro,ngro);
  ASSERT(ret == 0);


  /* PN : FIN Creation de la famille 0 */

  map < int, int > mapNoeud;
  typedef pair<set<int>::iterator, bool> IsFamily;
  int nbFamillesNoeud;


  i = 0;
  set<int> FamilySet;
  nbFamillesNoeud = 0;
  int verifienbnoeuds = 0;
  med_int* rien=0;

  SMDS_MeshNodesIterator itNodes(myMesh);
  for (;itNodes.More();itNodes.Next())
    {
      const Handle(SMDS_MeshElement)& elem = itNodes.Value();
      const Handle(SMDS_MeshNode)& node = myMesh->GetNode(1,elem);

      if(mdim==3){
	coo[i*3]=node->X();
	coo[i*3+1]=node->Y();
	coo[i*3+2]=node->Z();
      } else if(mdim==2) {
	if ( dimX ) {
	  coo[i*2]=node->Y();
	  coo[i*2+1]=node->Z();
	} 
	if ( dimY ) {
	  coo[i*2]=node->X();
	  coo[i*2+1]=node->Z();
	} 
	if ( dimZ ) {
	  coo[i*2]=node->X();
	  coo[i*2+1]=node->Y();
	} 
      } else {
	if ( dimX ) {
	  coo[i*2]=node->Y();
	  coo[i*2+1]=node->Z();
	} 
	if ( dimY ) {
	  coo[i*2]=node->X();
	  coo[i*2+1]=node->Z();
	} 
	if ( dimZ ) {
	  coo[i*2]=node->X();
	  coo[i*2+1]=node->Y();
	} 
      }
      mapNoeud[node->GetID()] = i+1;

      // renvoie 0 pour les noeuds internes du volume
      int numfamille=node->GetPosition()->GetShapeId();
      nufano[i]=numfamille;

      //SCRUTE(i);
      //SCRUTE(nufano[i]);
      //SCRUTE(coo[i*3]);
      //SCRUTE(coo[i*3+1]);
      //SCRUTE(coo[i*3+2]);
      if ( nufano[i] != 0 )
	{
	  IsFamily deja = FamilySet.insert(nufano[i]); // insert if new, or gives existant
	  if (deja.second)                       	      // actually inserted
	    {
	      char famille[MED_TAILLE_NOM+1];
	      sprintf(famille,"F%d",nufano[i]);
	      //	      CreateFamily(strdup(nommaa),strdup(famille),nufano[i],attvalabs++);
	      attvalabs++;
	      CreateFamily(strdup(nommaa),strdup(famille),nufano[i],numfamille);
	      //MESSAGE("---famille-noeud--- "<<nbFamillesNoeud<<" "<<nufano[i]);
	      nbFamillesNoeud++;
	    }
	}

      i++;
      verifienbnoeuds ++;
    }
  ret = MEDnoeudsEcr(myFileId,nommaa,mdim,coo,mode_coo,MED_CART,
		     nomcoo,unicoo,nomnoe,MED_FAUX,rien,MED_FAUX,
		     nufano,nnoe,MED_REMP);
  ASSERT(ret == 0);
  MESSAGE("--- Creation de " << verifienbnoeuds << " noeuds");
  ASSERT(verifienbnoeuds == nnoe);
  MESSAGE("--- Creation de " << nbFamillesNoeud << " familles de noeuds");

  /* liberation memoire */
  free(coo);
  free(numnoe);
  free(nufano);

  /****************************************************************************
   *                       ECRITURE DES ELEMENTS                              *
   ****************************************************************************/
  MESSAGE("(**************************)");
  MESSAGE("(* ELEMENTS DU MAILLAGE : *)");
  MESSAGE("(**************************)");
  
  /* Ecriture des connectivites, noms, numeros des mailles */

  if (ret == 0)
    {
      int nbFamillesElts =0;
      Handle(SMESHDS_Mesh) mySMESHDSMesh = Handle(SMESHDS_Mesh)::DownCast(myMesh);
      TopTools_IndexedMapOfShape myIndexToShape;
      TopExp::MapShapes(mySMESHDSMesh->ShapeToMesh(),myIndexToShape);

      map < int, int > mapFamille;

      if (besoinfamilledemaille == 1)
	{
	  int t;
	  for (t =1; t <= myIndexToShape.Extent(); t++)
	    {
	      const TopoDS_Shape S = myIndexToShape(t);
	      if (mySMESHDSMesh->HasMeshElements(S))
		{
		  //MESSAGE ("********* Traitement de la Famille "<<-t);
		
		  Handle(SMESHDS_SubMesh) SM = mySMESHDSMesh->MeshElements(S);
		  const TColStd_ListOfInteger& indElt = SM->GetIDElements();
		  TColStd_ListIteratorOfListOfInteger ite(indElt);
		
		  bool plein=false;
		  for (; ite.More(); ite.Next())
		    {
		      int eltId = ite.Value();
		      mapFamille [eltId] = - t;
		      plein=true;
		    }
		  if (plein)
		    {
		      nbFamillesElts++;
		      char famille[MED_TAILLE_NOM+1];
		      sprintf(famille,"E%d",t);
		      CreateFamily(strdup(nommaa),strdup(famille),-t,attvalabs++);
		    }
		}
	    }
	}

      int indice=1;
      for (i=0;i<MED_NBR_GEOMETRIE_MAILLE;i++)
	{
	  if (nmailles[i] > 0 && ret == 0)
	    {
	      MESSAGE ( " Start "<<typmai[i]);
	
	      /* dimension de la maille */
	      edim = typmai[i] / 100;
	      nsup = 0;
	      if (mdim  == 2 || mdim == 3)
		if (edim == 1) nsup = 1;
	      if (mdim == 3)
		if (edim == 2) nsup = 1;
	      //SCRUTE(nsup);
	
	      taille = nsup+typmai[i]%100;
	      //SCRUTE(taille);
	
	      /* allocation memoire */
	      connectivite = (med_int*)malloc(sizeof(med_int)* taille*nmailles[i]);
	      nomele = (char*)malloc(sizeof(char)*MED_TAILLE_PNOM* nmailles[i]+1);
	      numele = (med_int*)malloc(sizeof(med_int)* nmailles[i]);
	      nufael = (med_int*)malloc(sizeof(med_int)* nmailles[i]);
	      nomele="";
	      nbNodes=typmai[i]%100;
	
	      for (j=0;j<nmailles[i];j++)
		{
		  myId = elem_Id[i][j];
		  const Handle(SMDS_MeshElement)& elem = myMesh->FindElement(myId);
		  //*(numele+j) = myId;
		  *(numele+j) = indice ++;
		
		  for (k=0; k<nbNodes; k++)
		    {
		      *(connectivite+j*taille+k) = mapNoeud[elem->GetConnection(k+1)];
		      //SCRUTE(*(connectivite+j*taille+k));
		    }
		  if (nsup) *(connectivite+j*taille + nbNodes) = 0;
		
		  if (besoinfamilledemaille == 1)
		    {
		      if (mapFamille.find(myId)!=mapFamille.end())
			{
			  nufael[j]=mapFamille[myId];
			}
		      else
			{
			  nufael[j]=0;
			}
		    }
		  else
		    {
		      nufael[j]=0;
		    }
		
		  //SCRUTE(myId);
		  //SCRUTE(j);
		  //SCRUTE(nufael[j]);
		}
	
	      /* ecriture des données */
		
		    med_int* rien=0;
		    ret=MEDelementsEcr( myFileId,nommaa,mdim,connectivite,mode_coo,nomele,MED_FAUX,numele,
					MED_VRAI,nufael,nmailles[i],MED_MAILLE,typmai[i],typ_con,MED_REMP);
		    ASSERT(ret == 0);
		    //SCRUTE(ret);
		
		    if (ret < 0) MESSAGE(">> ERREUR : ecriture des mailles \n");
		
		    /* liberation memoire */
		    free(connectivite);
		    free(numele);
		    free(nufael);
		    MESSAGE ( " End "<<typmai[i]);
	    }
      };
	  MESSAGE("--- Creation de " << nbFamillesElts << " familles d elements");

    }

      /****************************************************************************
      *                      FERMETURE DU FICHIER                                *
      ****************************************************************************/

      ret = MEDfermer(myFileId);

      if (ret != 0)
	fprintf(stderr,">> ERREUR : erreur a la fermeture du fichier %s\n",file2Read);
      MESSAGE ("fichier ferme");


}

void DriverMED_W_SMESHDS_Mesh::CreateFamily(char* nommaa, char* famille, int i,med_int k) {

  med_int ngro=0;
  med_int natt;

  natt=1;
  char attdes[MED_TAILLE_DESC+1];
  char gro[MED_TAILLE_LNOM+1];
  char fam2[MED_TAILLE_LNOM+1];

  strcpy(attdes,"");
  strcpy(gro,"");
  strcpy(fam2,famille);

  med_int * attide=new med_int[1];
  med_int * attval=new med_int[1];
  attide[0] = k;
  attval[0] = k;


  //MESSAGE("-------- Creation de la Famille : "<< famille << "numero " << i << " --------------");
  med_int ret = MEDfamCr(myFileId, nommaa, fam2, i, attide, attval, attdes, natt, gro, ngro);
  ASSERT(ret==0);
  delete [] attide;
  delete [] attval;

}