vtf-logo

VisualGridData.h

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 // Copyright (C) 2002 Ralf Deiterding
00004 // Brandenburgische Universitaet Cottbus
00005 //
00006 // Copyright (C) 2003-2007 California Institute of Technology
00007 // Ralf Deiterding, ralf@amroc.net
00008 
00009 #ifndef AMROC_VISUALGRIDDATA_H
00010 #define AMROC_VISUALGRIDDATA_H
00011 
00019 #include "GridData3.h" 
00020 #include "Evaluator.h"
00021 
00022 #include <iostream>
00023 #include <vector>
00024 #include <cfloat>
00025 #include <cmath>
00026 
00027 template <class VectorType> class VisualGrid;
00028 
00035 template <class VectorType>
00036 class VisualGridData {
00037   typedef typename VectorType::InternalDataType DataType;
00038   typedef VisualGridData<VectorType> self;
00039 
00040 public:
00041   typedef VisualGrid<VectorType> grid_type;
00042   typedef Evaluator<VectorType> evaluator_type;
00043   typedef GridData<DataType,3> data_block_type;
00044   typedef GridData<VectorType,3> vector_block_type;
00045 
00046   typedef std::vector<self *> gblock_list_type;
00047   typedef typename gblock_list_type::iterator gblock_list_iterator;
00048   typedef std::vector<float> point_list_type;
00049   typedef typename point_list_type::iterator point_list_iterator;
00050 
00051   enum direction {x_dir, y_dir, z_dir};
00052 
00053 public:
00054   VisualGridData(grid_type& gr, data_block_type& db, BBox &bb, int comp, evaluator_type& ev,
00055                  int lev, int proc, float* pdx, float* pgeom) : 
00056     _Grid(gr), _Evaluator(ev), _Level(lev), _Processor(proc), dx(pdx), geom(pgeom) {
00057 
00058     _data_block = new vector_block_type(bb);    
00059     BeginFastIndex3(source, db.bbox(), db.data(), DataType);
00060     BeginFastIndex3(target, _data_block->bbox(), _data_block->data(),
00061                     VectorType);
00062     for_3 (i, j, k, bb, bb.stepsize())
00063       FastIndex3(target,i,j,k)(comp) = FastIndex3(source,i,j,k); 
00064     end_for
00065     EndFastIndex3(source);
00066     EndFastIndex3(target);
00067     DimUsed = gr.Dimensions();
00068   }
00069 
00070   VisualGridData(grid_type& gr, vector_block_type& db, BBox &bb, evaluator_type& ev,
00071                  int lev, int proc, float* pdx, float* pgeom) : 
00072     _Grid(gr), _Evaluator(ev), _Level(lev), _Processor(proc), dx(pdx), geom(pgeom) {
00073 
00074     _data_block = new vector_block_type(bb);
00075     BeginFastIndex3(source, db.bbox(), db.data(), VectorType);
00076     BeginFastIndex3(target, _data_block->bbox(), _data_block->data(),
00077                     VectorType);
00078     for_3 (i, j, k, bb, bb.stepsize())
00079       FastIndex3(target,i,j,k) = FastIndex3(source,i,j,k); 
00080     end_for
00081     EndFastIndex3(source);
00082     EndFastIndex3(target);
00083     DimUsed = gr.Dimensions();
00084   }
00085 
00086   virtual ~VisualGridData() 
00087   { delete _data_block; }
00088 
00089 public:
00090   virtual void SetBlocks(int blocks[],int& offset,int& nOffset) {
00091     std::cerr << ".";
00092     _nodeOffset = nOffset;
00093     blocks[offset]   = DB().extents(0)+DimUsed(0);
00094     blocks[offset+1] = DB().extents(1)+DimUsed(1);
00095     blocks[offset+2] = DB().extents(2)+DimUsed(2);
00096     nOffset += blocks[offset]*blocks[offset+1]*blocks[offset+2];   
00097     offset += 3;
00098   }
00099 
00100   virtual void SetGridCells(float xyz[][3],int& offset){
00101     std::cerr << ".";
00102     int i, j, k,orii,orij,orik,ssi,ssj,ssk,ijk,idx;
00103     int exti,extj,extk;
00104     exti=DB().extents(0);
00105     extj=DB().extents(1);
00106     extk=DB().extents(2);
00107     orii=DB().lower(0);
00108     orij=DB().lower(1);
00109     orik=DB().lower(2);
00110     ssi=DB().stepsize(0);
00111     ssj=DB().stepsize(1);
00112     ssk=DB().stepsize(2);
00113     
00114     idx = 0;
00115     for (k = 0; k < extk; k++){
00116       for (j = 0; j < extj; j++){
00117         for (i = 0; i < exti; i++){
00118           
00119           ijk = k*extj*exti + j*exti + i;
00120           
00121           xyz[ijk+offset][0] = (orii+((float)i+0.5)*ssi)*dx[0]+geom[0];
00122           xyz[ijk+offset][1] = (orij+((float)j+0.5)*ssj)*dx[1]+geom[2];
00123           xyz[ijk+offset][2] = (orik+((float)k+0.5)*ssk)*dx[2]+geom[4];
00124           idx++;
00125         }
00126       }
00127     }
00128     offset+=idx;
00129   }
00130   virtual void SetGridNodes(float xyz[][3],int& offset) {
00131     std::cerr << ".";
00132     int i, j, k,orii,orij,orik,ssi,ssj,ssk,ijk,idx;
00133     int exti,extj,extk;
00134     exti=DB().extents(0)+DimUsed(0);
00135     extj=DB().extents(1)+DimUsed(1);
00136     extk=DB().extents(2)+DimUsed(2);
00137     orii=DB().lower(0);
00138     orij=DB().lower(1);
00139     orik=DB().lower(2);
00140     ssi=DB().stepsize(0);
00141     ssj=DB().stepsize(1);
00142     ssk=DB().stepsize(2);
00143     
00144     idx = 0;
00145     for (k = 0; k < extk; k++){
00146       for (j = 0; j < extj; j++){
00147         for (i = 0; i < exti; i++){
00148           
00149           ijk = k*extj*exti + j*exti + i;
00150           
00151           xyz[ijk+offset][0] = (float)(orii+i*ssi)*dx[0]+geom[0];
00152           xyz[ijk+offset][1] = (float)(orij+j*ssj)*dx[1]+geom[2];
00153           xyz[ijk+offset][2] = (float)(orik+k*ssk)*dx[2]+geom[4];
00154           idx++;
00155         }
00156       }
00157     }
00158     offset+=idx;
00159   }
00160   
00161   virtual void SetGridConns(int con[][8],int& offset) {
00162     std::cerr << ".";
00163     int i,j,k,ijk,idx;
00164     int exti,extj,extk;
00165     exti=DB().extents(0)+DimUsed(0);
00166     extj=DB().extents(1)+DimUsed(1);
00167     extk=DB().extents(2)+DimUsed(2);
00168     
00169     idx = 0;
00170     for (k = 0; k < extk-DimUsed(2); k++){
00171       for (j = 0; j < extj-DimUsed(1); j++){
00172         for (i = 0; i < exti-DimUsed(0); i++){
00173           
00174           ijk = k*(extj-DimUsed(1))*(exti-DimUsed(0)) +
00175                 j*(exti-DimUsed(0)) + i;
00176           
00177           con[ijk+offset][0] =  k*extj*exti +
00178                                 j*exti +
00179                                 i + _nodeOffset;
00180           con[ijk+offset][1] =  k*extj*exti +
00181                                 j*exti +
00182                                 (i+DimUsed(0)) + _nodeOffset;
00183           con[ijk+offset][2] =  k*extj*exti +
00184                                 (j+DimUsed(1))*exti +
00185                                 i + _nodeOffset; 
00186           con[ijk+offset][3] =  k*extj*exti +
00187                                 (j+DimUsed(1))*exti +
00188                                 (i+DimUsed(0)) + _nodeOffset; 
00189           con[ijk+offset][4] = (k+DimUsed(2))*extj*exti +
00190                                j*exti +
00191                                i + _nodeOffset;
00192           con[ijk+offset][5] = (k+DimUsed(2))*extj*exti +
00193                                j*exti +
00194                                (i+DimUsed(0)) + _nodeOffset; 
00195           con[ijk+offset][6] = (k+DimUsed(2))*extj*exti +
00196                                (j+DimUsed(1))*exti +
00197                                i + _nodeOffset; 
00198           con[ijk+offset][7] = (k+DimUsed(2))*extj*exti +
00199                                (j+DimUsed(1))*exti +
00200                                (i+DimUsed(0)) + _nodeOffset; 
00201           
00202           idx++;
00203         }
00204       }
00205     }
00206     offset+=idx;
00207   }
00208   
00209   virtual void SetScalCells(int *key,float v[],int& offset,bool* CompRead) {
00210     std::cerr << ".";
00211     if (*key==1 || *key==2) {
00212       int idx=0;
00213       for_3 (i,j,k,DB().bbox(),DB().bbox().stepsize())
00214         if (*key==1) {
00215           v[idx+offset]= Processor() + 1.0;
00216           TheEvaluator().mvals[0] = Min(TheEvaluator().mvals[0],v[idx+offset]);
00217           TheEvaluator().mvals[1] = Max(TheEvaluator().mvals[1],v[idx+offset]);
00218         }
00219         else {
00220           v[idx+offset]= Level() + 1.0;   
00221           TheEvaluator().mvals[2] = Min(TheEvaluator().mvals[2],v[idx+offset]);
00222           TheEvaluator().mvals[3] = Max(TheEvaluator().mvals[3],v[idx+offset]);
00223         }
00224         idx++;        
00225       end_for
00226       offset+=idx;
00227     }
00228     else
00229       TheEvaluator().SetScalCells(*key,v,offset,DB(),dx,CompRead);
00230   }  
00231 
00232   virtual void SetScalNodes(int *key,float v[],int& offset,bool* CompRead) {
00233     std::cerr << ".";
00234     if (*key==1) {
00235       int idx=0;
00236       BBox lbbox(grow(DB().bbox(),DimUsed));
00237       data_block_type dbwork(lbbox);
00238       dbwork = -1.0;
00239       BeginFastIndex3(tgt, lbbox, dbwork.data(), DataType);
00240       
00241       for (gblock_list_iterator bit=_Grid.BlockList().begin();
00242            bit!=_Grid.BlockList().end(); bit++) 
00243         if ((*bit)->Processor() == Processor()) {
00244           BBox bbcur(growupper((*bit)->DB().bbox(),1));
00245           bbcur.setstepsize(lbbox.stepsize());
00246           bbcur.growupper(-1);
00247           BBox bbwork(lbbox * bbcur);
00248           if (bbwork.empty())
00249             continue;
00250           
00251           bbwork.lower() = (bbwork.lower()/bbwork.stepsize())*
00252                            bbwork.stepsize();
00253           bbwork.upper() = (bbwork.upper()/bbwork.stepsize())*
00254                            bbwork.stepsize();
00255           
00256           for_3 (i, j, k, bbwork, bbwork.stepsize())
00257             FastIndex3(tgt,i,j,k) = Processor();
00258           end_for
00259         } 
00260         
00261       int StepsizeMax = 1;
00262       for (int l=1; l<_Grid.NLevels(); l++)
00263         StepsizeMax *= _Grid.RefineFactor(l-1); 
00264 
00265       BBox bbox(DB().bbox());
00266       Coords step = DB().bbox().stepsize();
00267       Coords dimstep = DimUsed*step;
00268       bbox.upper() += dimstep;
00269       for_3 (i, j, k, bbox, step)
00270         v[idx+offset] = Processor() + 1.0;
00271       
00272         TheEvaluator().mvals[0] = Min(TheEvaluator().mvals[0], v[idx+offset]);
00273         TheEvaluator().mvals[1] = Max(TheEvaluator().mvals[1], v[idx+offset]);
00274         idx++;        
00275         end_for 
00276 
00277       EndFastIndex3(tgt); 
00278       offset+=idx;
00279     }
00280 
00281     else if (*key==2) {
00282       int idx=0;
00283       BBox lbbox(grow(DB().bbox(),DimUsed));
00284       data_block_type dbwork(lbbox);
00285       dbwork = -1.0;
00286       BeginFastIndex3(tgt, lbbox, dbwork.data(), DataType);
00287       
00288       for (gblock_list_iterator bit=_Grid.BlockList().begin();
00289            bit!=_Grid.BlockList().end(); bit++) 
00290         if ((*bit)->Level() == Level()) {
00291           BBox bbcur(growupper((*bit)->DB().bbox(),1));
00292           bbcur.setstepsize(lbbox.stepsize());
00293           bbcur.growupper(-1);
00294           BBox bbwork(lbbox * bbcur);
00295           if (bbwork.empty())
00296             continue;
00297           
00298           bbwork.lower() = (bbwork.lower()/bbwork.stepsize())*
00299                            bbwork.stepsize();
00300           bbwork.upper() = (bbwork.upper()/bbwork.stepsize())*
00301                            bbwork.stepsize();
00302           
00303           for_3 (i, j, k, bbwork, bbwork.stepsize())
00304             FastIndex3(tgt,i,j,k) = Level();
00305           end_for
00306         } 
00307         
00308       int StepsizeMax = 1;
00309       for (int l=1; l<_Grid.NLevels(); l++)
00310         StepsizeMax *= _Grid.RefineFactor(l-1); 
00311 
00312       BBox bbox(DB().bbox());
00313       Coords step = DB().bbox().stepsize();
00314       Coords dimstep = DimUsed*step;
00315       bbox.upper() += dimstep;
00316       for_3 (i, j, k, bbox, step)
00317         v[idx+offset] = Level() + 1.0;
00318       
00319         TheEvaluator().mvals[2] = Min(TheEvaluator().mvals[2], v[idx+offset]);
00320         TheEvaluator().mvals[3] = Max(TheEvaluator().mvals[3], v[idx+offset]);
00321         idx++;        
00322         end_for 
00323 
00324       EndFastIndex3(tgt); 
00325       offset+=idx;
00326     }
00327 
00328     else {
00329       BBox lbbox(grow(DB().bbox(),DimUsed));
00330       lbbox *= _Grid.GlobalBBox();
00331       lbbox.lower() = (lbbox.lower()/lbbox.stepsize())*lbbox.stepsize();
00332       lbbox.upper() = (lbbox.upper()/lbbox.stepsize())*lbbox.stepsize();
00333       vector_block_type dbwork(lbbox);
00334       dbwork = FLT_MAX;
00335       BeginFastIndex3(tgt, lbbox, dbwork.data(), VectorType);
00336         
00337       for (gblock_list_iterator bit=_Grid.BlockList().begin();
00338            bit!=_Grid.BlockList().end(); bit++) {
00339         BBox bbcur(growupper((*bit)->DB().bbox(),1));
00340         bbcur.setstepsize(lbbox.stepsize());
00341         bbcur.growupper(-1);
00342         BBox bbwork(lbbox * bbcur);
00343         if (bbwork.empty())
00344           continue;
00345         
00346         bbwork.lower() = (bbwork.lower()/bbwork.stepsize())*bbwork.stepsize();
00347         bbwork.upper() = (bbwork.upper()/bbwork.stepsize())*bbwork.stepsize();
00348 
00349         Coords step_src = (*bit)->DB().bbox().stepsize();
00350         Coords step_tgt = bbwork.stepsize();
00351         if (step_tgt(0) == step_src(0)) 
00352           dbwork.copy((*bit)->DB(), bbwork);
00353         
00354         else if (step_tgt(0) > step_src(0)) {
00355           BeginFastIndex3(src, (*bit)->DB().bbox(), (*bit)->DB().data(),
00356                           VectorType);
00357           for_3 (i, j, k, bbwork, bbwork.stepsize())
00358             BBox AddOver(3,i,j,k,
00359                          i+step_tgt(0)-step_src(0),
00360                          j+step_tgt(1)-step_src(1),
00361                          k+step_tgt(2)-step_src(2),
00362                          step_src(0),step_src(1),step_src(2));
00363             AddOver *= (*bit)->DB().bbox();
00364             float c=0;
00365             FastIndex3(tgt,i,j,k) = 0.0;
00366             for_3 (l, m, n, AddOver, step_src)
00367               FastIndex3(tgt,i,j,k) += FastIndex3(src,l,m,n);
00368               c += 1.0;
00369             end_for
00370             if (c>0) FastIndex3(tgt,i,j,k) /= c;
00371           end_for
00372           EndFastIndex3(src);
00373         }
00374 
00375         else {
00376           BeginFastIndex3(src, (*bit)->DB().bbox(), (*bit)->DB().data(),
00377                           VectorType);
00378           for_3 (i, j, k, bbwork, bbwork.stepsize())
00379             FastIndex3(tgt,i,j,k) = FastIndex3(src,i-i%step_src(0),
00380                                                j-j%step_src(1),
00381                                                k-k%step_src(2));
00382           end_for
00383           EndFastIndex3(src);
00384         }
00385       } 
00386       EndFastIndex3(tgt); 
00387         
00388       BBox bbox(DB().bbox());
00389       Coords dimstep = DimUsed*DB().bbox().stepsize();
00390       bbox.upper() += dimstep;
00391       TheEvaluator().SetScalNodes(*key,v,offset,dbwork,bbox,dx,CompRead);
00392     }
00393   }
00394 
00395   virtual void SetPointsCells(int *key,point_list_type &p,float v[],bool* CompRead) {
00396     register int d;
00397     float lb[3], ub[3];
00398     for (d=0; d<3; d++) {
00399       lb[d]=DB().lower(d)*dx[d]+geom[2*d];
00400       ub[d]=(DB().extents(d)+DimUsed(d)-1)*DB().stepsize(d)*dx[d]+lb[d];
00401     }
00402     std::vector<unsigned int> plist;
00403     for (register unsigned int i=0;i<p.size()/3;i++) {
00404       for (d=0; d<3; d++) 
00405         if (DimUsed(d))
00406           if (p[3*i+d]<lb[d] || p[3*i+d]>ub[d]) break;
00407       if (d==3) plist.push_back(i); 
00408     }
00409 
00410     if (!plist.empty()) {
00411       float *tmpdata = new float[NCells()];
00412       int offset=0;
00413       SetScalCells(key,tmpdata,offset,CompRead);
00414       for (register unsigned int i=0;i<plist.size();i++) {
00415         int idx[3];
00416         for (d=0; d<3; d++) {
00417           if (DimUsed(d))
00418             idx[d]=int((p[3*plist[i]+d]-lb[d])/(DB().stepsize(d)*dx[d]));
00419           else idx[d]=0;
00420           assert (idx[d]>=0 && idx[d]<DB().extents(d));
00421           // Do not test points anymore after data has been assigned once
00422           p[3*plist[i]+d]=-1.e37; 
00423         }
00424         v[plist[i]]=tmpdata[idx[0]+DB().extents(0)*(idx[1]+DB().extents(1)*idx[2])];
00425       }
00426       delete [] tmpdata;
00427     }
00428   }
00429 
00430   virtual void SetPointsNodes(int *key,point_list_type &p,float v[],bool* CompRead) {
00431     register int d;
00432     float lb[3], ub[3];
00433     for (d=0; d<3; d++) {
00434       lb[d]=DB().lower(d)*dx[d]+geom[2*d];
00435       ub[d]=(DB().extents(d)+DimUsed(d)-1)*DB().stepsize(d)*dx[d]+lb[d];
00436     }
00437     std::vector<unsigned int> plist;
00438     for (register unsigned int i=0;i<p.size()/3;i++) {
00439       for (d=0; d<3; d++) 
00440         if (DimUsed(d))
00441           if (p[3*i+d]<lb[d] || p[3*i+d]>ub[d]) break;
00442       if (d==3) plist.push_back(i); 
00443     }
00444 
00445     if (!plist.empty()) {
00446       float *tmpdata = new float[NNodes()];
00447       int offset=0;
00448       SetScalNodes(key,tmpdata,offset,CompRead);
00449       for (register unsigned int i=0;i<plist.size();i++) {
00450         int idx[3], idp[3], ext[3];
00451         float fac[3];
00452         for (d=0; d<3; d++) {
00453           if (DimUsed(d)) {
00454             idx[d]=int((p[3*plist[i]+d]-lb[d])/(DB().stepsize(d)*dx[d]));
00455             fac[d]=(p[3*plist[i]+d]-lb[d]-idx[d]*DB().stepsize(d)*dx[d])/
00456               (DB().stepsize(d)*dx[d]);
00457             idp[d]=idx[d]+1;
00458           }
00459           else {
00460             idx[d]=0; fac[d]=0.; idp[d]=0;
00461           }
00462           ext[d]=DB().extents(d)+DimUsed(d);
00463           assert (idx[d]>=0 && idp[d]<ext[d]);
00464           // Do not test points anymore after data has been assigned once
00465           p[3*plist[i]+d]=-1.e37; 
00466         }
00467         
00468         float &d0=tmpdata[idx[0]+ext[0]*(idx[1]+ext[1]*idx[2])];
00469         float &d1=tmpdata[idp[0]+ext[0]*(idx[1]+ext[1]*idx[2])];
00470         float &d2=tmpdata[idx[0]+ext[0]*(idp[1]+ext[1]*idx[2])];
00471         float &d3=tmpdata[idp[0]+ext[0]*(idp[1]+ext[1]*idx[2])];
00472         float &d4=tmpdata[idx[0]+ext[0]*(idx[1]+ext[1]*idp[2])];
00473         float &d5=tmpdata[idp[0]+ext[0]*(idx[1]+ext[1]*idp[2])];
00474         float &d6=tmpdata[idx[0]+ext[0]*(idp[1]+ext[1]*idp[2])];
00475         float &d7=tmpdata[idp[0]+ext[0]*(idp[1]+ext[1]*idp[2])];
00476         
00477         v[plist[i]]=(1.0-fac[2])*((1.-fac[1])*((1.-fac[0])*d0+fac[0]*d1)+
00478                                       fac[1] *((1.-fac[0])*d2+fac[0]*d3))+
00479                          fac[2] *((1.-fac[1])*((1.-fac[0])*d4+fac[0]*d5)+
00480                                       fac[1] *((1.-fac[0])*d6+fac[0]*d7));
00481       }
00482       delete [] tmpdata;
00483     }
00484   }
00485 
00486   virtual void VectCells(int *key,float v[][3],int& offset,bool* CompRead) {
00487     std::cerr << ".";
00488     int idx=0;  
00489     BeginFastIndex3(dat, DB().bbox(), DB().data(), VectorType);
00490     for_3 (i,j,k,DB().bbox(),DB().bbox().stepsize())
00491       v[idx+offset][0]= FastIndex3(dat,i,j,k)(1);
00492       v[idx+offset][1]= FastIndex3(dat,i,j,k)(2);
00493       v[idx+offset][2]= FastIndex3(dat,i,j,k)(3);           
00494       idx++;          
00495     end_for
00496     offset+=idx;
00497     EndFastIndex3(dat); 
00498   }
00499 
00500   virtual void VectNodes(int *key,float v[][3],int& offset,bool* CompRead) {
00501     std::cerr << ".";
00502     int idx=0;
00503 
00504     BBox lbbox(grow(DB().bbox(),DimUsed));
00505     lbbox *= _Grid.GlobalBBox();
00506     lbbox.lower() = (lbbox.lower()/lbbox.stepsize())*lbbox.stepsize();
00507     lbbox.upper() = (lbbox.upper()/lbbox.stepsize())*lbbox.stepsize();
00508     vector_block_type dbwork(lbbox);
00509     dbwork = FLT_MAX;
00510     BeginFastIndex3(tgt, lbbox, dbwork.data(), VectorType);
00511     
00512     for (gblock_list_iterator bit=_Grid.BlockList().begin();
00513          bit!=_Grid.BlockList().end(); bit++) {
00514       BBox bbcur(growupper((*bit)->DB().bbox(),1));
00515       bbcur.setstepsize(lbbox.stepsize());
00516       bbcur.growupper(-1);
00517       BBox bbwork(lbbox * bbcur);
00518       if (bbwork.empty())
00519         continue;
00520 
00521       bbwork.lower() = (bbwork.lower()/bbwork.stepsize())*bbwork.stepsize();
00522       bbwork.upper() = (bbwork.upper()/bbwork.stepsize())*bbwork.stepsize();
00523       
00524       Coords step_src = (*bit)->DB().bbox().stepsize();
00525       Coords step_tgt = bbwork.stepsize();
00526       if (step_tgt(0) == step_src(0)) 
00527         dbwork.copy((*bit)->DB(), bbwork);
00528 
00529       else if (step_tgt(0) > step_src(0)) {
00530         BeginFastIndex3(src, (*bit)->DB().bbox(), (*bit)->DB().data(),
00531                         VectorType);
00532         for_3 (i, j, k, bbwork, bbwork.stepsize())
00533           BBox AddOver(3,i,j,k,
00534                        i+step_tgt(0)-step_src(0),
00535                        j+step_tgt(1)-step_src(1),
00536                        k+step_tgt(2)-step_src(2),
00537                        step_src(0),step_src(1),step_src(2));
00538           AddOver *= (*bit)->DB().bbox();
00539           float c=0;
00540           FastIndex3(tgt,i,j,k) = 0.0;
00541           for_3 (l, m, n, AddOver, step_src)
00542             FastIndex3(tgt,i,j,k) += FastIndex3(src,l,m,n);
00543             c += 1.0;
00544           end_for
00545           if (c>0) FastIndex3(tgt,i,j,k) /= c;
00546         end_for
00547         EndFastIndex3(src);
00548       }
00549 
00550       else {
00551         BeginFastIndex3(src, (*bit)->DB().bbox(), (*bit)->DB().data(),
00552                         VectorType);
00553         for_3 (i, j, k, bbwork, bbwork.stepsize())
00554           FastIndex3(tgt,i,j,k) = FastIndex3(src,i-i%step_src(0),
00555                                              j-j%step_src(1),
00556                                              k-k%step_src(2));
00557         end_for
00558         EndFastIndex3(src);
00559       }
00560     } 
00561       
00562     BBox bbox(DB().bbox());
00563     Coords dimstep = DimUsed*DB().bbox().stepsize();
00564     bbox.upper() += dimstep;
00565       
00566     for_3 (i, j, k, bbox, dimstep)
00567       BBox AddOver(3,i-dimstep(0),j-dimstep(1),k-dimstep(2),
00568                    i,j,k,dimstep(0),dimstep(1),dimstep(2));
00569       AddOver *= lbbox;
00570           
00571       float c=0;
00572       v[idx+offset][0]=0;
00573       v[idx+offset][1]=0;
00574       v[idx+offset][2]=0;
00575       for_3 (l, m, n, AddOver, dimstep)
00576         v[idx+offset][0] += FastIndex3(tgt,l,m,n)(1);
00577         v[idx+offset][1] += FastIndex3(tgt,l,m,n)(2);
00578         v[idx+offset][2] += FastIndex3(tgt,l,m,n)(3);
00579         c += 1.0;
00580       end_for
00581       if (c>0) {
00582         v[idx+offset][0] /= c;
00583         v[idx+offset][1] /= c;
00584         v[idx+offset][2] /= c;
00585       }
00586       idx++;          
00587     end_for 
00588       
00589     EndFastIndex3(tgt); 
00590     offset+=idx;
00591   }
00592 
00593   virtual vector_block_type& DB() { return *_data_block; }
00594   virtual int NNodes() { return (!DB().bbox().empty() ? 
00595      ((DB().extents(0)+DimUsed(0))*
00596       (DB().extents(1)+DimUsed(1))*
00597       (DB().extents(2)+DimUsed(2))) : 0); }
00598   virtual int NCells() { return DB().size(); }
00599   virtual const int& Level() const { return _Level; } 
00600   virtual int Level() { return _Level; } 
00601   virtual const int& Processor() const { return _Processor; } 
00602   virtual int Processor() { return _Processor; } 
00603   virtual evaluator_type& TheEvaluator() { return _Evaluator; }
00604 
00605 protected:
00606   grid_type& _Grid;
00607   vector_block_type* _data_block;
00608   evaluator_type& _Evaluator;
00609   int _Level;
00610   int _Processor;
00611   float* dx;
00612   float* geom;
00613   int _nodeOffset;
00614   Coords DimUsed;
00615 };
00616 
00617 #endif

Generated on Fri Aug 24 13:00:31 2007 for AMROC's HDF Tools - by  doxygen 1.4.7