1 /*  This file is part of MED.
2  *
3  *  COPYRIGHT (C) 1999 - 2019  EDF R&D, CEA/DEN
4  *  MED is free software: you can redistribute it and/or modify
5  *  it under the terms of the GNU Lesser General Public License as published by
6  *  the Free Software Foundation, either version 3 of the License, or
7  *  (at your option) any later version.
8  *
9  *  MED is distributed in the hope that it will be useful,
10  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  *  GNU Lesser General Public License for more details.
13  *
14  *  You should have received a copy of the GNU Lesser General Public License
15  *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
16  */
17 
18 
19 #include <med.h>
20 #include <med_config.h>
21 #include <med_outils.h>
22 #include <hdf5.h>
23 
24 
25 #ifdef TRACEPFL
26 #define tracepfl(x) x
27 #else
28 #define tracepfl(x)
29 #endif
30 
31 /*
32  * - Nom de la fonction : _MEDdatasetNumLire
33  * - Description : lecture d'un dataset tableau numerique
34  * - Parametres :
35  *     - pere (IN)     : l'ID de l'objet HDF pere ou placer l'attribut
36  *     - nom  (IN)     : le nom du dataset
37  *     - type (IN)     : type numerique MED
38  *     - interlace (IN) : Choix du type d'entrelacement demandé par l'appelant { MED_FULL_INTERLACE(x1,y1,z1,x2,...)) , MED_NO_INTERLACE(x1,x2,y1,y2,z1,z2) }
39  *       - nbdim   (IN) : Dimension des éléments
40  *       - fixdim  (IN) : MED_ALL ou n° de la dimension a enregistrer à partir de 1..oo
41  *     - psize     (IN) : Taille du profil à utiliser, MED_NOPF si pas de profil
42  *       - pflmod  (IN) : Indique comment lire les informations en mémoire { MED_COMPACT, MED_GLOBAL }.
43  *       - pflcmp  (IN) : MED_PFL_NON_COMPACT = fichier profil non compacté, MED_PFL_COMPACT = fichier profil compacté
44  *       - pfltab  (IN) : Tableau contenant les n° déléments à traiter (1....oo)
45  *       - ngauss  (IN) : Nombre de points de GAUSS par élément
46  *     - nbelem (IN)    : Taille du tableau de valeurs
47  *                        (référence tous les élements, cette taille  prend en compte le nombre de pts de gauss et la dimension )
48  *                        (n'est utile que pour la lecture en mode MED_GLOBAL d'un champ avec profil compacté : pflcmp = MED_PFL_COMPACT)
49  *     - val  (OUT)    : valeurs du tableau
50  * - Resultat : 0 en cas de succes, -1 sinon
51  *  Equivalent à l'ancienne routine si .....,MED_NO_INTERLACE,1,MED_ALL,MED_NOPF,0,1 (peu importe),....
52  */
53 
54 
_MEDdatasetNumLire232(int dummy,...)55 void  _MEDdatasetNumLire232(int dummy,...) {
56 
57   va_list params;
58 
59   med_idt              pere;
60   char                *nom;
61   med_type_champ       type;
62   med_mode_switch interlace;
63   med_size        nbdim,fixdim,psize,*pfltab;
64   med_mode_profil pflmod;
65   med_int         ngauss,nbelem;
66   med_lecture_profil pflcmp;
67   unsigned char   *val;
68   med_err         *fret;
69 
70   med_idt    dataset, dataspace = 0, memspace = 0;
71   med_size   start_mem[1],start_data[1],*pflmem=0,*pfldsk=0;
72   med_size   stride[1],count[1],pcount[1],size[1],pflsize[1];
73   med_err    ret;
74   med_idt    i,j,index;
75   med_idt    type_hdf;
76   hid_t      datatype;
77   size_t     typesize;
78   int        dim, firstdim, dimutil, lastdim;
79   med_size   sizencmp[1];
80 
81   va_start(params,dummy);
82   pere      = va_arg(params,med_idt);
83   nom       = va_arg(params,char *);
84   type      = va_arg(params,med_type_champ);
85   interlace = va_arg(params,med_mode_switch);
86   nbdim     = va_arg(params,med_size);
87   fixdim    = va_arg(params,med_size);
88   psize     = va_arg(params,med_size);
89   pflmod    = va_arg(params,med_mode_profil);
90   pflcmp    = va_arg(params,med_lecture_profil);
91   pfltab    = va_arg(params,med_size *);
92   ngauss    = va_arg(params,med_int);
93   nbelem    = va_arg(params,med_int);
94   val       = va_arg(params,  unsigned char *);
95   fret      = va_arg(params,  med_err *);
96 
97 
98   /* Verify fixdim is between [0, nbdim] ( 0 is MED_ALL ) */
99   if ( fixdim > nbdim  )
100     goto Fail;
101 
102   switch(type)
103     {
104     case MED_FLOAT64 :
105       type_hdf = H5T_NATIVE_DOUBLE;
106       break;
107 
108     case MED_INT32 :
109       type_hdf = H5T_NATIVE_INT;
110       break;
111 
112     case MED_INT64 :
113       type_hdf = H5T_NATIVE_LONG;
114       break;
115 
116     default :
117       goto Fail;
118     }
119 
120   /* Ouverture du Dataset à lire */
121   if ((dataset = H5Dopen(pere,nom)) < 0)
122     goto Fail;
123 
124   /* Interrogation de la taille du dataset */
125   if ( (datatype  = H5Dget_type(dataset )) < 0) goto Fail;
126   if ( (typesize  = H5Tget_size(datatype)) == 0) goto Fail;
127 #if 1
128   {
129   hid_t   space;
130   hsize_t   sizespace[H5S_MAX_RANK];
131   hsize_t   maxsizespace[H5S_MAX_RANK];
132   int       ndims;
133 
134   space = H5Dget_space(dataset);
135   ndims = H5Sget_simple_extent_dims(space, sizespace, maxsizespace);
136   H5Sclose(space);
137   size[0] = sizespace[0];
138   }
139 #else
140   size[0] = H5Dget_storage_size(dataset) / typesize;
141 #endif
142   if ( H5Tclose(datatype) < 0) goto Fail;
143 
144 
145   /* verification (facultative) de la taille du dataset */
146   if ( psize != MED_NOPF) {
147     if ( pflcmp == MED_PFL_COMPACT ) {
148       if (size[0] != psize * ngauss * nbdim)
149         goto Fail;
150     }
151   }
152 
153   /* Create dataspace */
154   if ((dataspace = H5Screate_simple(1,size,NULL)) < 0)
155     goto Fail;
156 
157   switch(interlace) {
158     case MED_FULL_INTERLACE :
159 
160       /*Initialisation des indices de boucle du traitement de l'entrelacement en fonction de la dimension fixee*/
161       if ( fixdim != MED_ALL) {
162 	    firstdim = fixdim-1;
163 	    lastdim  = fixdim;
164 	    dimutil  = 1;
165 	  }
166       else	{
167 	    firstdim = 0;
168 	    lastdim = nbdim;
169 	    dimutil  = nbdim;
170 	  }
171 
172       count [0] = (*size)/(nbdim);
173 
174 
175       /*rem: Pas de vérification de l'assertion (*size)=n*nbdim */
176       if ( psize == MED_NOPF ) {
177 
178         tracepfl(printf("%s branche 1 : %lld\n", __FILE__, fixdim));
179 
180 	    /* Creation d'un data space mémoire de dimension 1, de longeur size, et de longeur maxi size */
181 	    if ( (memspace = H5Screate_simple (1, size, NULL)) <0)
182 	      goto Fail;
183 
184 	    stride[0] = nbdim;
185 
186 	    for (dim=firstdim; dim < lastdim; dim++) {
187 
188 	      start_mem[0] = dim;
189 	      if ( (ret = H5Sselect_hyperslab (memspace, H5S_SELECT_SET, start_mem, stride,
190 					   count, NULL)) <0)
191 	        goto Fail;
192 
193 	      start_data[0] = dim*count[0];
194 	      if ( (ret = H5Sselect_hyperslab (dataspace, H5S_SELECT_SET, start_data, NULL,
195 					   count, NULL)) <0)
196 	        goto Fail;
197 
198 	      if ((ret = H5Dread(dataset,type_hdf,memspace,dataspace,
199 			     H5P_DEFAULT, val)) < 0)
200 	        goto Fail;
201 	    }
202 
203       }
204       else {
205 
206 	    pflsize [0] = psize*ngauss*nbdim;
207 	    pcount  [0] = psize*ngauss*dimutil;
208 	    pflmem     = (med_size *) malloc (sizeof(med_size)*pcount[0]);
209 	    pfldsk     = (med_size *) malloc (sizeof(med_size)*pcount[0]);
210 
211 	    switch(pflmod) { /* switch pflmod pour FULL_INTERLACE*/
212 	      case MED_GLOBAL :
213 
214             tracepfl(printf("%s branche 2 : %lld - %d\n", __FILE__, fixdim, pflcmp));
215 
216 	        /* Creation d'un data space mémoire de dimension 1, de longeur size, et de longeur maxi size */
217             if (pflcmp == MED_PFL_COMPACT)
218               sizencmp[0] = nbelem * nbdim;
219             else
220               sizencmp[0] = size[0];
221             if ( (memspace = H5Screate_simple (1, sizencmp, NULL)) <0)
222               goto Fail;
223 
224 	        for (dim=firstdim; dim < lastdim; dim++) {
225 
226 	          for (i=0; i < psize; i++)              /* i balaye les élements du profil */
227 	            for (j=0; j < ngauss; j++) {
228 		          index = i*ngauss+j + (dim-firstdim)*(psize*ngauss);
229 		          pflmem[index] = (pfltab[i]-1)*ngauss*nbdim + j*nbdim+dim;
230 		          if (pflcmp == MED_PFL_COMPACT)
231 		            pfldsk[index] = dim*count[0] + i*ngauss+j;
232 		          else
233   		            pfldsk[index] = dim*count[0] + (pfltab[i]-1)*ngauss+j;
234 	            }
235 	        }
236 
237 	        if ( (ret = H5Sselect_elements(memspace ,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pflmem ) ) <0)
238 	          goto Fail;
239 
240 	        if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
241 	          goto Fail;
242 
243 	        break;
244 
245 	      case MED_COMPACT :
246 
247             tracepfl(printf("%s branche 3 : %lld - %d\n", __FILE__, fixdim, pflcmp));
248 
249 	        /* Creation d'un data space mémoire de dimension 1, de la longeur du profil          */
250 	        /* La dimension utilisée est ici nbdim, même pour un profil compact on suppose       */
251 	        /*  que l'utilisateur a toutes les coordonées stockées, même si il en demande qu'une */
252 
253 	        if ( (memspace = H5Screate_simple (1, pflsize, NULL)) <0)
254 	          goto Fail;
255 
256 	        for (dim=firstdim; dim < lastdim; dim++) {
257 
258 	          for (i=0; i < psize; i++)              /* i balaye les élements du profil */
259 	            for (j=0; j < ngauss; j++) {
260 		          index = i*ngauss+j + (dim-firstdim)*(psize*ngauss);
261 		          pflmem[index] = i*ngauss*nbdim + j*nbdim+dim;
262 		          if (pflcmp == MED_PFL_COMPACT)
263 		            pfldsk[index] = dim*count[0] + i*ngauss+j;
264 		          else
265   		            pfldsk[index] = dim*count[0] + (pfltab[i]-1)*ngauss+j;
266 	            }
267 	        }
268 
269 	        if ( (ret = H5Sselect_elements(memspace ,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pflmem ) ) <0)
270 	          goto Fail;
271 
272 	        if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
273 	          goto Fail;
274 
275 	        break;
276 
277 	      default :
278 	        goto Fail;
279 	    }
280 
281 	    if ((ret = H5Dread(dataset,type_hdf,memspace,dataspace,H5P_DEFAULT, val)) < 0)
282 	      goto Fail;
283 
284 	    free(pflmem);
285 	    free(pfldsk);
286       }
287 
288       break;
289 
290     case MED_NO_INTERLACE :
291 
292       /*Initialisation des indices de boucle du traitement de l'entrelacement en fonction de la dimension fixee*/
293 
294       count[0] = (*size)/nbdim;
295 
296       if ( psize == MED_NOPF ) {
297 
298         tracepfl(printf("%s branche 4 : %d\n", __FILE__, fixdim));
299 
300 	    if ( fixdim != MED_ALL)
301 	      start_data[0] = (fixdim-1)*count[0];
302 	    else {
303 	      count[0] = *size;
304 	      start_data[0] =  0;
305 	    };
306 
307 	    if ( (ret = H5Sselect_hyperslab (dataspace, H5S_SELECT_SET, start_data, NULL,
308 					 count, NULL)) <0)
309 	      goto Fail;
310 
311 	    if ((ret = H5Dread(dataset,type_hdf,dataspace,dataspace,
312 			   H5P_DEFAULT, val)) < 0)
313 	      goto Fail;
314 
315       }
316       else {
317 
318 	    if ( fixdim != MED_ALL) {
319 	      firstdim = fixdim-1;
320 	      lastdim  = fixdim;
321 	      dimutil  = 1;
322 	    }
323         else	{
324 	      firstdim = 0;
325 	      lastdim  = nbdim;
326 	      dimutil  = nbdim;
327 	    }
328 
329 	    pflsize [0] = psize*ngauss*nbdim;   /* taille du memspace (toutes les composantes) */
330   	    pcount  [0] = psize*ngauss*dimutil; /* taille des indexes de selection, pflmem et pfldsk  */
331 	    /*nom pas très coherent avec count !!! A revoir */
332 	    pfldsk      = (med_size *) malloc(sizeof(med_size)*pcount[0]);
333 
334 	    switch(pflmod) { /*switch plfmod pour NO_INTERLACE */
335 	      case MED_GLOBAL :
336 
337             tracepfl(printf("%s branche 5 : %lld - %d\n", __FILE__, fixdim, pflcmp));
338 
339 	        if ( pflcmp == MED_PFL_COMPACT ) {
340 	          pflmem     = (med_size *) malloc (sizeof(med_size)*pcount[0]);
341 	          sizencmp[0] = nbelem * nbdim;
342               if ( (memspace = H5Screate_simple (1, sizencmp, NULL)) <0)
343   	            goto Fail;
344   	        }
345 
346 	        for (dim=firstdim; dim < lastdim; dim++) {
347 
348  	          for (i=0; i < psize; i++)              /* i balaye le nbre d'élements du profil                */
349 		        for (j=0; j < ngauss; j++) {
350 		          index = i*ngauss+j + (dim-firstdim)*(psize*ngauss);
351 		          if ( pflcmp == MED_PFL_COMPACT ) {
352 		            pflmem[index] = dim*nbelem+(pfltab[i]-1)*ngauss+j;
353 		            pfldsk[index] = dim*count[0] + i*ngauss+j;
354 		          }
355 		          else {
356   		            pfldsk[index] = dim*count[0] + (pfltab[i]-1)*ngauss+j;
357   		          }
358 		        }
359 	        }
360 
361 	        if ( pflcmp == MED_PFL_COMPACT ) {
362 	          if ( (ret = H5Sselect_elements(memspace ,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pflmem ) ) <0)
363 	            goto Fail;
364 
365 	          if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET,pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
366 	            goto Fail;
367 
368 	          if ((ret = H5Dread(dataset,type_hdf,memspace,dataspace,H5P_DEFAULT, val)) < 0)
369 	            goto Fail;
370 
371 	          free(pflmem);
372             }
373             else {
374 
375 	          if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET,pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
376 	            goto Fail;
377 
378 	          if ((ret = H5Dread(dataset,type_hdf,dataspace,dataspace,H5P_DEFAULT, val)) < 0)
379 	            goto Fail;
380 	        }
381 
382 	        break;
383 
384 	      case MED_COMPACT :
385 
386             tracepfl(printf("%s branche 6 : %lld - %d\n", __FILE__, fixdim, pflcmp));
387 
388 	        if ( pflcmp == MED_PFL_NON_COMPACT ) {
389               /* Creation d'un data space mémoire de dimension 1, de la longeur du profil          */
390 	          /* La dimension utilisée est ici nbdim, même pour un profil compact on suppose       */
391 	          /*  que l'utilisateur a toutes les coordonées stockées, même si il en demande qu'une */
392 
393 	          if ( (memspace = H5Screate_simple (1, pflsize, NULL)) <0)
394 	            goto Fail;
395 
396 	          /*La taille de pflmem n'est pas forcément égale à celle de memspace, ex : choix de lire 1 composante*/
397 	          pflmem     = (med_size *) malloc (sizeof(med_size)*pcount[0]);
398 	        }
399 
400 	        /* Le profil COMPACT est contigüe, mais il est possible que l'on selectionne uniquemenent une dimension*/
401 
402 	        index = 0;
403 	        for (dim=firstdim; dim < lastdim; dim++) {
404 	          for (i=0; i < psize; i++)              /* i balaye le nbre d'élements du profil */
405 		        for (j=0; j < ngauss; j++) {
406 		          /*FORMULATION 1ere : index = i*ngauss+j + (dim-firstdim)*(psize*ngauss);*/
407 		          /*FORMULATION 2sd  : index = ( (dim-firstdim)*psize + i )*ngauss + j;   */
408 		          /*FORMULATION 1ere : pflmem[index] = dim*(psize*ngauss) + i*ngauss+j;*/
409 		          if ( pflcmp == MED_PFL_NON_COMPACT ) {
410                     pflmem[index] = ( (dim*psize) + i )*ngauss + j;
411 		            pfldsk[index] = dim*count[0]  + (pfltab[i]-1)*ngauss+j;
412 		          }
413                   else {
414                     pfldsk[index] = dim*count[0] + i*ngauss+j;
415                   }
416 		          index++;
417 		        }
418 	        }
419 
420 	        if ( pflcmp == MED_PFL_NON_COMPACT ) {
421               if ( (ret = H5Sselect_elements(memspace ,H5S_SELECT_SET, pcount[0], HDF5_SELECT_BUG pflmem ) ) <0)
422 	            goto Fail;
423 
424 	          if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET,pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
425 	            goto Fail;
426 
427 	          if ((ret = H5Dread(dataset,type_hdf,memspace,dataspace,H5P_DEFAULT, val)) < 0)
428 	            goto Fail;
429 
430 	          free(pflmem);
431 	        }
432             else {
433 	          if ( (ret = H5Sselect_elements(dataspace,H5S_SELECT_SET,pcount[0], HDF5_SELECT_BUG pfldsk ) ) <0)
434 	            goto Fail;
435 
436 	          if ((ret = H5Dread(dataset,type_hdf,dataspace,dataspace,H5P_DEFAULT, val)) < 0)
437 	            goto Fail;
438             }
439 
440 	        break;
441 
442 	      default :
443 	        goto Fail;
444 
445 	    }
446 
447 	    free(pfldsk);
448 
449       };
450 
451       break;
452 
453     default :
454       goto Fail;
455   }
456 
457 
458 
459   if (memspace)
460     if ((ret = H5Sclose(memspace)) < 0)
461       goto Fail;
462 
463   if ((ret = H5Sclose(dataspace)) < 0)
464     goto Fail;
465 
466   if ((ret = H5Dclose(dataset)) < 0)
467     goto Fail;
468 
469  Success:
470   va_end(params);
471   *fret=0;
472   return;
473 
474  Fail:
475   va_end(params);
476   *fret = -1;
477   return;
478 }
479