00001 #include <stdio.h>
00002 #include <stdlib.h>
00003 #include <string.h>
00004 #include <hdf5.h>
00005 #include "H5Part.h"
00006
00007 static unsigned h5part_debug=0;
00008
00009
00010 herr_t H5PartIOcounter(hid_t group_id,
00011 const char *member_name,
00012 void *operator_data);
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033 H5PartFile *H5PartOpenFileParallel(const char *filename,unsigned flags
00034 #ifdef PARALLEL_IO
00035 , MPI_Comm comm
00036 #endif
00037 ){
00038 H5PartFile *f=(H5PartFile*)malloc(sizeof(H5PartFile));
00039 #ifdef PARALLEL_IO
00040 MPI_Info info = MPI_INFO_NULL;
00041 #endif
00042 f->xfer_prop = f->create_prop = f->access_prop=H5P_DEFAULT;
00043 #ifdef PARALLEL_IO
00044 MPI_Comm_size(comm,&f->nprocs);
00045 MPI_Comm_rank(comm,&f->myproc);
00046 f->pnparticles=(long long*)malloc(sizeof(long long)*f->nprocs);
00047
00048 f->access_prop = H5Pcreate(H5P_FILE_ACCESS);
00049 if(H5Pset_fapl_mpio(f->access_prop,comm,info)<0){
00050 fprintf(stderr,"Total failure trying to setup mpi-io for access\n");
00051 exit(0);
00052 }
00053
00054
00055 f->create_prop = H5P_DEFAULT;
00056
00057
00058
00059 f->xfer_prop = H5Pcreate(H5P_DATASET_XFER);
00060 H5Pset_dxpl_mpio(f->xfer_prop,H5FD_MPIO_COLLECTIVE);
00061 f->comm = comm;
00062 #endif
00063 if(flags==H5PART_READ){
00064 f->file=H5Fopen(filename,H5F_ACC_RDONLY,f->access_prop);
00065
00066
00067 }
00068 else if(flags==H5PART_WRITE){
00069 f->file=H5Fcreate(filename,H5F_ACC_TRUNC,f->create_prop,f->access_prop);
00070 }
00071 else {
00072 f->file=-1;
00073 fprintf(stderr,"H5Open: Invalid file access type %d\n",flags);
00074 }
00075 if(f->file<=0) {
00076 free(f);
00077 fprintf(stderr,"H5Open: File open failed for file [%s]\n",
00078 filename);
00079 exit(0);
00080 return 0;
00081 }
00082 #ifdef PARALLEL_IO
00083 else {
00084 if(h5part_debug) fprintf(stderr,"Proc[%u] Opened file %s val=%d\n",
00085 f->myproc,filename,f->file);
00086 }
00087 #endif
00088 f->mode=flags;
00089 f->maxstep=0;
00090 f->timegroup=0;
00091 f->timestep=0;
00092 f->shape=0;
00093 f->diskshape=H5S_ALL;
00094 f->memshape=H5S_ALL;
00095 f->viewstart=-1;
00096 f->viewend=-1;
00097 return f;
00098 }
00099
00100 H5PartFile *H5PartOpenFile(const char *filename,unsigned flags){
00101 H5PartFile *f=(H5PartFile*)malloc(sizeof(H5PartFile));
00102
00103
00104 f->xfer_prop = f->create_prop = f->access_prop=H5P_DEFAULT;
00105 #ifdef PARALLEL_IO
00106 f->pnparticles=0;
00107 f->comm = MPI_COMM_WORLD;
00108 f->nprocs=1;
00109 f->myproc=0;
00110 #endif
00111 if(flags==H5PART_READ){
00112
00113
00114 f->file = H5Fopen(filename,H5F_ACC_RDONLY,H5P_DEFAULT);
00115 if(f->file==-1){
00116 fprintf(stderr,"H5Open: File %s not found\n",filename);
00117 }
00118 }
00119 else if(flags==H5PART_WRITE){
00120 f->file=H5Fcreate(filename,H5F_ACC_TRUNC,H5P_DEFAULT,H5P_DEFAULT);
00121 }
00122 else {
00123 f->file=-1;
00124 fprintf(stderr,"H5Open: Invalid file access type %d\n",flags);
00125 }
00126 if(f->file<0) {
00127 free(f);
00128 fprintf(stderr,"H5Open: File open failed for file [%s]\n",
00129 filename);
00130 return 0;
00131 }
00132 f->mode=flags;
00133 f->maxstep=0;
00134 f->timegroup=0;
00135 f->timestep=0;
00136 f->shape=0;
00137 f->diskshape=H5S_ALL;
00138 f->memshape=H5S_ALL;
00139 f->viewstart=-1;
00140 f->viewend=-1;
00141 return f;
00142 }
00143 int H5PartFileIsValid(H5PartFile *f){
00144 if(f->file > 0) return 1;
00145 else return 0;
00146 }
00147
00148
00149
00150 void H5PartCloseFile(H5PartFile *f){
00151 if(f->shape>0) H5Sclose(f->shape); f->shape=0;
00152 if(f->timegroup>0) H5Gclose(f->timegroup); f->timegroup=0;
00153 if(f->diskshape!=H5S_ALL) H5Sclose(f->diskshape);
00154 #ifdef PARALLEL_IO
00155 if(f->pnparticles) free(f->pnparticles);
00156 #endif
00157 if(f->xfer_prop!=H5P_DEFAULT) {
00158 H5Pclose(f->xfer_prop); f->xfer_prop=H5P_DEFAULT;
00159 }
00160 if(f->access_prop!=H5P_DEFAULT) {
00161 H5Pclose(f->access_prop); f->access_prop=H5P_DEFAULT;
00162 }
00163 if(f->create_prop!=H5P_DEFAULT) {
00164 H5Pclose(f->create_prop); f->create_prop=H5P_DEFAULT;
00165 }
00166 if(f->file) H5Fclose(f->file); f->file=0;
00167 free(f);
00168 }
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179 void H5PartSetNumParticles(H5PartFile *f,long long nparticles){
00180 hssize_t start[1];
00181 hsize_t stride[1],count[1],total;
00182 register int i,r;
00183
00184 #ifndef PARALLEL_IO
00185
00186
00187
00188
00189 if(f->nparticles==nparticles) return;
00190 #endif
00191 if(f->diskshape!=H5S_ALL) H5Sclose(f->diskshape);
00192 if(f->memshape!=H5S_ALL) H5Sclose(f->memshape);
00193 f->memshape=f->diskshape = H5S_ALL;
00194 if(f->shape) H5Sclose(f->shape);
00195 f->nparticles=(hsize_t)nparticles;
00196 #ifndef PARALLEL_IO
00197 f->shape=H5Screate_simple(1,
00198 &(f->nparticles),
00199 NULL);
00200 #else
00201
00202
00203
00204
00205
00206
00207
00208
00209
00210 MPI_Allgather(&nparticles,1,MPI_LONG_LONG,f->pnparticles,1,MPI_LONG_LONG,f->comm);
00211 if(h5part_debug) {
00212 if(f->myproc==0){
00213 puts("AllGather: Particle offsets:\n");
00214 for(i=0;i<f->nprocs;i++)
00215 printf("\tnp=%u\n",f->pnparticles[i]);
00216 }
00217 }
00218
00219 stride[0]=1;
00220 start[0]=0;
00221 for(i=0;i<f->myproc;i++) start[0]+=f->pnparticles[i];
00222 total=0;
00223 for(i=0;i<f->nprocs;i++) total+=f->pnparticles[i];
00224 f->shape = H5Screate_simple(1,&total,NULL);
00225 f->diskshape = H5Screate_simple(1,&total,NULL);
00226 f->memshape = H5Screate_simple(1,&(f->nparticles),NULL);
00227
00228 if(f->shape<0 || f->memshape<0 || f->diskshape<0){
00229 fprintf(stderr,"Abort: shape construction failed\n");
00230 exit(0);
00231 }
00232 count[0]=nparticles;
00233
00234 r=H5Sselect_hyperslab(f->diskshape,H5S_SELECT_SET,start,stride,count,NULL);
00235 if(r<0){
00236 fprintf(stderr,"Abort: Selection Failed!\n");
00237 exit(0);
00238 }
00239 if(f->timegroup<0){
00240 H5PartSetStep(f,0);
00241 }
00242 #endif
00243 }
00244
00245
00246
00247
00248
00249
00250
00251
00252
00253
00254 int H5PartWriteDataFloat64(H5PartFile *f,char *name,double *array){
00255 register int r;
00256 hid_t dataset;
00257 if(f->mode==H5PART_READ){
00258 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00259 return 0;
00260 }
00261
00262
00263 if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
00264 dataset = H5Dcreate(f->timegroup,name,H5T_NATIVE_DOUBLE,f->shape,H5P_DEFAULT);
00265 if(dataset<0)
00266 fprintf(stderr,"Dataset open failed for name=[%s] step=%u!\n",
00267 name,f->timestep);
00268 r=H5Dwrite(dataset,H5T_NATIVE_DOUBLE,f->memshape,f->diskshape,H5P_DEFAULT,array);
00269 H5Dclose(dataset);
00270 return r;
00271 }
00272 int H5PartWriteDataInt64(H5PartFile *f,char *name,long long *array){
00273 register int r;
00274 hid_t dataset;
00275
00276
00277 if(f->mode==H5PART_READ){
00278 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00279 return 0;
00280 }
00281 if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
00282 dataset = H5Dcreate(f->timegroup,name,H5T_NATIVE_INT64,f->shape,H5P_DEFAULT);
00283 if(dataset<0){
00284 fprintf(stderr,"Dataset open failed for name=[%s] step=%u!\n",
00285 name,f->timestep);
00286 exit(0);
00287 }
00288 r=H5Dwrite(dataset,H5T_NATIVE_INT64,f->memshape,f->diskshape,H5P_DEFAULT,array);
00289 if(r<0) {
00290 fprintf(stderr,"Attempt to write dataset failed!\n");
00291 exit(0);
00292 }
00293 H5Dclose(dataset);
00294 return r;
00295 }
00296 int H5PartWriteFileAttribString(H5PartFile *f,char *name,
00297 char *attrib){
00298 if(f->mode==H5PART_READ){
00299 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00300 return 0;
00301 }
00302 return H5PartWriteFileAttrib(f,name,H5T_NATIVE_CHAR,attrib,strlen(attrib)+1);
00303 }
00304 int H5PartWriteStepAttribString(H5PartFile *f,char *name,
00305 char *attrib){
00306 if(f->mode==H5PART_READ){
00307 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00308 return 0;
00309 }
00310 return H5PartWriteStepAttrib(f,name,H5T_NATIVE_CHAR,attrib,strlen(attrib)+1);
00311 }
00312 int H5PartWriteStepAttrib(H5PartFile *f,char *name,
00313 hid_t type,void *value,int nelem){
00314 register int r;
00315 hid_t attrib;
00316 hid_t space;
00317 hsize_t len;
00318 if(f->mode==H5PART_READ){
00319 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00320 return 0;
00321 }
00322
00323
00324 if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
00325 len = nelem;
00326 space = H5Screate_simple(1,&len,NULL);
00327 attrib = H5Acreate(f->timegroup,name,type,space,H5P_DEFAULT);
00328 r=H5Awrite(attrib,type,value);
00329 H5Aclose(attrib);
00330 H5Sclose(space);
00331 return r;
00332 }
00333
00334 int H5PartWriteAttrib(H5PartFile *f,char *name,
00335 hid_t type,void *value,int nelem){
00336 if(f->mode==H5PART_READ){
00337 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00338 return 0;
00339 }
00340 return H5PartWriteStepAttrib(f,name,type,value,nelem);
00341 }
00342
00343
00344 int H5PartWriteFileAttrib(H5PartFile *f,char *name,
00345 hid_t type,void *value,int nelem){
00346 register int r;
00347 hid_t attrib;
00348 hid_t space;
00349 hsize_t len;
00350 if(f->mode==H5PART_READ){
00351 fprintf(stderr,"H5PartWriteDataFloat64: Error! Attempting to write to read-only file\n");
00352 return 0;
00353 }
00354
00355
00356 if(f->file<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!! File is not open!\n");}
00357 len = nelem;
00358 space = H5Screate_simple(1,&len,NULL);
00359 attrib = H5Acreate(f->file,name,type,space,H5P_DEFAULT);
00360 r=H5Awrite(attrib,type,value);
00361 H5Aclose(attrib);
00362 H5Sclose(space);
00363 return r;
00364 }
00365 herr_t H5PartAttribcounter(hid_t group_id,
00366 const char *member_name,
00367 void *operator_data){
00368 int *count = (int*)operator_data;
00369 (*count)++;
00370 return 0;
00371 }
00372 int H5PartGetNumStepAttribs(H5PartFile *f){
00373 return H5Aget_num_attrs(f->timegroup);
00374 }
00375 int H5PartGetNumFileAttribs(H5PartFile *f){
00376
00377
00378 int nattribs;
00379 hid_t rootgroup=H5Gopen(f->file,"/");;
00380 nattribs=H5Aget_num_attrs(rootgroup);
00381
00382
00383 H5Gclose(rootgroup);
00384 return nattribs;
00385 }
00386
00387 hid_t H5PartNormType(hid_t type){
00388 H5T_class_t tclass = H5Tget_class(type);
00389 int size = H5Tget_size(type);
00390 switch(tclass){
00391 case H5T_INTEGER:
00392 if(size==8) return H5T_NATIVE_INT64;
00393 else if(size==1) return H5T_NATIVE_CHAR;
00394 else fprintf(stderr,"Unknown type %u. Cannot infer normalized type\n",type);
00395 break;
00396 case H5T_FLOAT:
00397 return H5T_NATIVE_DOUBLE;
00398 break;
00399 }
00400 return -1;
00401 }
00402
00403 void H5PartGetStepAttribInfo(H5PartFile *f,int idx,
00404 char *name,size_t maxname,hid_t *type,int *nelem){
00405 hid_t attrib=H5Aopen_idx(f->timegroup,idx);
00406 hid_t mytype=H5Aget_type(attrib);
00407 hid_t space = H5Aget_space(attrib);
00408 if(nelem)
00409 *nelem=H5Sget_simple_extent_npoints(space);
00410 if(name)
00411 H5Aget_name(attrib,maxname,name);
00412 if(type)
00413 *type=H5PartNormType(mytype);
00414
00415
00416
00417
00418 H5Sclose(space);
00419 H5Tclose(mytype);
00420 H5Aclose(attrib);
00421 }
00422 void H5PartGetFileAttribInfo(H5PartFile *f,int idx,
00423 char *name,size_t maxname,hid_t *type,int *nelem){
00424 hid_t attrib=H5Aopen_idx(f->file,idx);
00425 hid_t mytype=H5Aget_type(attrib);
00426 hid_t space = H5Aget_space(attrib);
00427 if(nelem)
00428 *nelem=H5Sget_simple_extent_npoints(space);
00429 if(name)
00430 H5Aget_name(attrib,maxname,name);
00431 if(type)
00432 *type=H5PartNormType(mytype);
00433
00434
00435
00436
00437 H5Sclose(space);
00438 H5Tclose(mytype);
00439 H5Aclose(attrib);
00440 }
00441
00442 int H5PartReadStepAttrib(H5PartFile *f,char *name,void *data){
00443
00444 hid_t attrib=H5Aopen_name(f->timegroup,name);
00445 if(attrib<=0) return -1;
00446 hid_t mytype=H5Aget_type(attrib);
00447 hid_t space = H5Aget_space(attrib);
00448 hsize_t nelem=H5Sget_simple_extent_npoints(space);
00449 hid_t type=H5PartNormType(mytype);
00450
00451
00452
00453
00454 H5Aread(attrib,type,data);
00455 H5Sclose(space);
00456 H5Tclose(mytype);
00457 H5Aclose(attrib);
00458 return 1;
00459 }
00460
00461 void H5PartReadAttrib(H5PartFile *f,char *name,void *data){
00462
00463 H5PartReadStepAttrib(f,name,data);
00464 }
00465
00466
00467
00468
00469 int H5PartReadFileAttrib(H5PartFile *f,char *name,void *data){
00470 hid_t attrib=H5Aopen_name(f->file,name);
00471 if(attrib<0) return -1;
00472 hid_t mytype=H5Aget_type(attrib);
00473 hid_t space = H5Aget_space(attrib);
00474 hsize_t nelem=H5Sget_simple_extent_npoints(space);
00475 hid_t type=H5PartNormType(mytype);
00476
00477
00478
00479
00480 H5Aread(attrib,type,data);
00481 H5Sclose(space);
00482 H5Tclose(mytype);
00483 H5Aclose(attrib);
00484 return 1;
00485 }
00486
00487
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499 void H5PartSetStep(H5PartFile *f,int step){
00500 char name[128];
00501 #ifdef PARALLEL_IO
00502 if(h5part_debug) printf("Proc[%u] SetStep to %u for file %u\n",f->myproc,step,f->file);
00503 #else
00504 if(h5part_debug) fprintf(stderr,"SetStep to %d for file\n",step);
00505 #endif
00506 f->timestep=step;
00507 if(f->timegroup>0) {
00508 H5Gclose(f->timegroup);
00509 }
00510 f->timegroup=-1;
00511 sprintf(name,"Particles#%u",f->timestep);
00512
00513
00514 if(f->mode==H5PART_READ) {
00515
00516 f->timegroup=H5Gopen(f->file,name);
00517 }
00518 if(f->timegroup<=0){
00519
00520
00521 if(f->mode == H5PART_READ){
00522 fprintf(stderr,"Error in H5PartSetStep() Step #%d does not exist!\n",step);
00523 }
00524 f->timegroup = H5Gcreate(f->file,name,0);
00525 if(f->timegroup<0) {
00526 #ifdef PARALLEL_IO
00527 fprintf(stderr,"p[%u] timegroup creation failed!\n",
00528 f->myproc);
00529 #endif
00530 return;
00531 }
00532 }
00533
00534
00535
00536
00537
00538 }
00539
00540
00541
00542
00543
00544
00545
00546
00547 herr_t H5PartIOcounter(hid_t group_id,
00548 const char *member_name,
00549 void *operator_data){
00550 int *count = (int*)operator_data;
00551 H5G_stat_t objinfo;
00552
00553 if(!strncmp(member_name,"Particles",9)) (*count)++;
00554 return 0;
00555 }
00556 herr_t H5PartDScounter(hid_t group_id,
00557 const char *member_name,
00558 void *operator_data){
00559 int *count = (int*)operator_data;
00560 H5G_stat_t objinfo;
00561
00562 if(H5Gget_objinfo(group_id,
00563 member_name,
00564 1 ,
00565 &objinfo)<0) {
00566 return 0;
00567 }
00568 if(objinfo.type==H5G_DATASET) (*count)++;
00569 return 0;
00570 }
00571
00572 typedef struct H5IO_getname_t {
00573 int index,count;
00574 char *name;
00575 }H5IO_getname_t;
00576
00577 herr_t H5IOgetname(hid_t group_id,
00578 const char *member_name,
00579 void *operator_data){
00580 H5IO_getname_t *getn = (H5IO_getname_t*)operator_data;
00581
00582 H5G_stat_t objinfo;
00583
00584 if(H5Gget_objinfo(group_id,
00585 member_name,
00586 1 ,
00587 &objinfo)<0) return 0;
00588
00589 if(objinfo.type!=H5G_DATASET)
00590 return 0;
00591 if(getn->index==getn->count){
00592 strcpy(getn->name,member_name);
00593 return 1;
00594 }
00595 getn->count++;
00596 return 0;
00597 }
00598
00599
00600
00601
00602
00603
00604
00605 int H5PartGetNumSteps(H5PartFile *f){
00606
00607 int count=0;
00608 int idx=0;
00609 while(H5Giterate(f->file,
00610 "/",
00611 &idx,
00612 H5PartIOcounter,
00613 &count)<0){}
00614 return count;
00615 }
00616 int H5PartGetNumDatasets(H5PartFile *f){
00617 int count=0;
00618 int idx=0;
00619 char name[128];
00620 sprintf(name,"Particles#%u",f->timestep);
00621 while(H5Giterate(f->file,
00622 name,
00623 &idx,
00624 H5PartDScounter,
00625 &count)<0){}
00626 return count;
00627 }
00628 int H5PartGetDatasetName(H5PartFile *f,
00629 int _index,
00630 char *name,
00631 int maxlen){
00632 H5IO_getname_t getn;
00633 int idx=_index;
00634 char gname[128];
00635 sprintf(gname,"Particles#%u",f->timestep);
00636 getn.index=_index; getn.name=name; getn.count=_index;
00637 while(H5Giterate(f->file,
00638 gname,
00639 &idx,
00640 H5IOgetname,
00641 &getn)<0){}
00642 return 1;
00643 }
00644
00645 hid_t H5PartGetDiskShape(H5PartFile *f,hid_t dataset){
00646 hid_t space = H5Dget_space(dataset);
00647 if(h5part_debug) fprintf(stderr,"H5PartGetDiskShape\n");
00648 if(H5PartHasView(f)){
00649 int r;
00650 hsize_t total, stride, count;
00651 hssize_t range[2];
00652 if(h5part_debug) fprintf(stderr,"\tSelection is available\n");
00653
00654 range[0]=f->viewstart;
00655 range[1]=f->viewend;
00656 count = range[1]-range[0]+1;
00657 stride=1;
00658
00659 if(f->diskshape>0)
00660 r=H5Sselect_hyperslab(f->diskshape,H5S_SELECT_SET,
00661 range,
00662 &stride,&count,NULL);
00663
00664 r=H5Sselect_hyperslab(space,H5S_SELECT_SET,
00665 range,&stride,&count,NULL);
00666 if(h5part_debug) fprintf(stderr,"\tSelection: range=%u:%u, npoints=%u s=%u\n",
00667 (int)range[0],(int)range[1],
00668 (int)H5Sget_simple_extent_npoints(space),
00669 (int)H5Sget_select_npoints(space));
00670 if(r<0){
00671 fprintf(stderr,"Abort: Selection Failed!\n");
00672 return space;
00673 }
00674 }
00675 else {
00676 if(h5part_debug) fprintf(stderr,"\tNo Selection\n");
00677 }
00678 return space;
00679 }
00680
00681 hid_t H5PartGetMemShape(H5PartFile *f,hid_t dataset){
00682 if(h5part_debug) fprintf(stderr,"H5PartGetMemShape\n");
00683 if(H5PartHasView(f)) {
00684 hsize_t len = f->viewend - f->viewstart+1;
00685 return H5Screate_simple(1,&len,NULL);
00686 }
00687 else {
00688 return H5S_ALL;
00689 }
00690 }
00691
00692
00693
00694
00695
00696
00697 long long H5PartGetNumParticles(H5PartFile *f){
00698 hid_t space,dataset;
00699 hsize_t nparticles;
00700 if(f->timegroup<0) {
00701 if(f->timestep<0)
00702 H5PartSetStep(f,0);
00703 else
00704 H5PartSetStep(f,f->timestep);
00705 }
00706 if(f->timegroup<=0){
00707 fprintf(stderr,"Damnit!!! tried to select a timestep %d\n", f->timestep);
00708 exit(0);
00709 }
00710
00711 dataset=H5Dopen(f->timegroup,"x");
00712
00713
00714 space = H5PartGetDiskShape(f,dataset);
00715 if(H5PartHasView(f)){
00716 if(h5part_debug) fprintf(stderr,"\tGetNumParticles::get_selection\n");
00717 nparticles=H5Sget_select_npoints(space);
00718 }
00719 else {
00720 if(h5part_debug) fprintf(stderr,"\tGetNumParticles::get_simple_extent from space\n");
00721 nparticles= H5Sget_simple_extent_npoints(space);
00722 }
00723 if(space!=H5S_ALL) H5Sclose(space);
00724 H5Dclose(dataset);
00725 return (long long)nparticles;
00726 }
00727
00728
00729
00730
00731
00732
00733
00734
00735
00736
00737
00738
00739
00740
00741
00742
00743 void H5PartSetView(H5PartFile *f,long long start,long long end){
00744 hsize_t total, stride, count;
00745 hssize_t range[2];
00746 int r;
00747 if(h5part_debug) fprintf(stderr,"SetView(%d,%d)\n",
00748 (int)start,(int)end);
00749 if(f->mode==H5PART_WRITE){
00750 fprintf(stderr,"H5PartSetView(): SetView does not make sense for write-only files. It is meant to be used for read-only files. (maybe later this will change)\n");
00751 return;
00752 }
00753
00754 f->viewstart=-1;
00755 f->viewend=-1;
00756 if(f->shape!=0){
00757 H5Sclose(f->shape);
00758 f->shape=0;
00759 }
00760 if(f->diskshape!=0 && f->diskshape!=H5S_ALL){
00761 H5Sclose(f->diskshape);
00762 f->diskshape=H5S_ALL;
00763 }
00764 f->diskshape = H5S_ALL;
00765 if(f->memshape!=0 && f->memshape!=H5S_ALL){
00766 H5Sclose(f->memshape);
00767 f->memshape=H5S_ALL;
00768 }
00769 if(start==-1 && end==-1) {
00770 if(h5part_debug) fprintf(stderr,"\tEarly Termination: Unsetting View\n");
00771 return;
00772 }
00773
00774
00775 total=H5PartGetNumParticles(f);
00776 if(h5part_debug) fprintf(stderr,"\tTotal nparticles=%u\n",(int)total);
00777 if(start==-1) start=0;
00778 if(end==-1) end=total-1;
00779
00780
00781
00782
00783
00784 range[0]=start;
00785 range[1]=end;
00786
00787 f->viewstart=range[0];
00788 f->viewend=range[1];
00789 f->nparticles=range[1]-range[0]+1;
00790 if(h5part_debug) fprintf(stderr,"\tRange is now %u:%u\n",(int)range[0],(int)range[1]);
00791
00792
00793
00794 f->shape = H5Screate_simple(1,&total,NULL);
00795
00796 f->diskshape= H5Screate_simple(1,&total,NULL);
00797 f->memshape = H5Screate_simple(1,&(f->nparticles),NULL);
00798
00799 if(f->shape<0 || f->memshape<0 || f->diskshape<0){
00800 fprintf(stderr,"Abort: shape construction failed\n");
00801 exit(0);
00802 }
00803 count = range[1]-range[0];
00804 if(h5part_debug) fprintf(stderr,"\tcount=%d\n",(int)count);
00805 stride=1;
00806
00807 r=H5Sselect_hyperslab(f->diskshape,H5S_SELECT_SET,range,&stride,&count,NULL);
00808 if(r<0){
00809 fprintf(stderr,"Abort: Selection Failed!\n");
00810 exit(0);
00811 }
00812
00813 }
00814
00815
00816
00817
00818
00819 int H5PartGetView(H5PartFile *f,long long *start,long long *end){
00820 long long range[2];
00821 range[0]=(f->viewstart>=0)?f->viewstart:0;
00822 range[1]=(f->viewend>=0)?f->viewend:H5PartGetNumParticles(f);
00823 if(start) {
00824 *start=range[0];
00825 }
00826 if(end) {
00827 *end=range[1];
00828 }
00829
00830 return range[1]-range[0];
00831 }
00832
00833
00834
00835
00836
00837
00838
00839
00840 void H5PartSetCanonicalView(H5PartFile *f){
00841
00842
00843
00844 if(f->mode != H5PART_READ){
00845 fprintf(stderr,"H5PartSetCanonicalView(): Views do not make sense for write-only files. It is meant to be used for read-only files. (maybe later this will change)\n");
00846 return;
00847 }
00848 H5PartSetView(f,-1,-1);
00849 #ifdef PARALLEL_IO
00850 if(f->timegroup<0){
00851 H5PartSetStep(f,0);
00852 }
00853
00854
00855
00856
00857
00858 if(H5PartReadStepAttrib(f,"pnparticles",f->pnparticles)<0){
00859
00860 int i;
00861 long long total=0, n = H5PartGetNumParticles(f);
00862 n/=f->nprocs;
00863 f->pnparticles[0]=n;
00864 total=n;
00865 for(i=1;i<f->nprocs;i++){
00866 f->pnparticles[i]=n;
00867 total+=n;
00868 }
00869 }
00870 {
00871 int i;
00872 long long total = 0, n = H5PartGetNumParticles(f);
00873
00874 for(i=0;i<f->myproc;i++){
00875 total+=f->pnparticles[i];
00876 }
00877 H5PartSetView(f,total,total+f->pnparticles[f->myproc]-1);
00878 }
00879 #endif
00880
00881
00882 }
00883
00884
00885
00886
00887
00888
00889
00890
00891 int H5PartReadDataFloat64(H5PartFile *f,char *name,double *array){
00892 hid_t space,memspace,dataset,datatype;
00893 if(!f->timegroup) H5PartSetStep(f,f->timestep);
00894 dataset=H5Dopen(f->timegroup,name);
00895 space = H5PartGetDiskShape(f,dataset);
00896 memspace = H5PartGetMemShape(f,dataset);
00897
00898 H5Dread(dataset,
00899 H5T_NATIVE_DOUBLE,
00900
00901 memspace,
00902 space,
00903 H5P_DEFAULT,
00904 array);
00905 if(space!=H5S_ALL) H5Sclose(space);
00906 if(memspace!=H5S_ALL) H5Sclose(memspace);
00907 H5Dclose(dataset);
00908 return 1;
00909 }
00910
00911 int H5PartReadDataInt64(H5PartFile *f,char *name,long long *array){
00912 hid_t space,memspace,dataset,datatype;
00913 if(!f->timegroup) H5PartSetStep(f,f->timestep);
00914 dataset=H5Dopen(f->timegroup,name);
00915 space = H5PartGetDiskShape(f,dataset);
00916 memspace = H5PartGetMemShape(f,dataset);
00917
00918 H5Dread(dataset,
00919 H5T_NATIVE_INT64,
00920
00921 memspace,
00922 space,
00923 H5P_DEFAULT,
00924 array);
00925 if(space!=H5S_ALL) H5Sclose(space);
00926 if(memspace!=H5S_ALL) H5Sclose(memspace);
00927 H5Dclose(dataset);
00928 return 1;
00929 }
00930
00931
00932
00933
00934
00935
00936
00937
00938
00939
00940 int H5PartReadParticleStep(H5PartFile *f,
00941 int step,
00942 double *x,double *y,double *z,
00943 double *px,double *py,double *pz,
00944 long long *id){
00945 H5PartSetStep(f,step);
00946
00947 H5PartReadDataFloat64(f,"x",x);
00948 H5PartReadDataFloat64(f,"y",y);
00949 H5PartReadDataFloat64(f,"z",z);
00950 H5PartReadDataFloat64(f,"px",px);
00951 H5PartReadDataFloat64(f,"py",py);
00952 H5PartReadDataFloat64(f,"pz",pz);
00953 H5PartReadDataInt64(f,"id",id);
00954 return 1;
00955 }
00956