/[escript]/trunk/escript/src/Data.cpp
ViewVC logotype

Diff of /trunk/escript/src/Data.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2245 by jfenwick, Thu Jan 8 23:55:40 2009 UTC revision 2246 by jfenwick, Thu Feb 5 06:04:55 2009 UTC
# Line 541  Data::copyWithMask(const Data& other, Line 541  Data::copyWithMask(const Data& other,
541    {    {
542      throw DataException("Error - Unknown DataAbstract passed to copyWithMask.");      throw DataException("Error - Unknown DataAbstract passed to copyWithMask.");
543    }    }
544      unsigned int selfrank=getDataPointRank();
545      unsigned int otherrank=other2.getDataPointRank();
546      unsigned int maskrank=mask2.getDataPointRank();
547      if ((selfrank==0) && (otherrank>0 || maskrank>0))
548      {
549        // to get here we must be copying from a large object into a scalar
550        // I am not allowing this.
551        // If you are calling copyWithMask then you are considering keeping some existing values
552        // and so I'm going to assume that you don't want your data objects getting a new shape.
553        throw DataException("Attempt to copyWithMask into a scalar from an object or mask with rank>0.");
554      }
555    // Now we iterate over the elements    // Now we iterate over the elements
556    DataVector& self=getReadyPtr()->getVector();    DataVector& self=getReadyPtr()->getVector();
557    const DataVector& ovec=other2.getReadyPtr()->getVector();    const DataVector& ovec=other2.getReadyPtr()->getVector();
558    const DataVector& mvec=mask2.getReadyPtr()->getVector();    const DataVector& mvec=mask2.getReadyPtr()->getVector();
559    if ((self.size()!=ovec.size()) || (self.size()!=mvec.size()))  
560      if ((selfrank>0) && (otherrank==0) &&(maskrank==0))
561      {
562        // Not allowing this combination.
563        // it is not clear what the rank of the target should be.
564        // Should it be filled with the scalar (rank stays the same);
565        // or should the target object be reshaped to be a scalar as well.
566        throw DataException("Attempt to copyWithMask from scalar mask and data into non-scalar target.");
567      }
568      if ((selfrank>0) && (otherrank>0) &&(maskrank==0))
569      {
570        if (mvec[0]>0)      // copy whole object if scalar is >0
571        {
572            copy(other);
573        }
574        return;
575      }
576      if (isTagged())       // so all objects involved will also be tagged
577      {
578        // note the !
579        if (!((getDataPointShape()==mask2.getDataPointShape()) &&
580            ((other2.getDataPointShape()==mask2.getDataPointShape()) || (otherrank==0))))
581        {
582            throw DataException("copyWithMask, shape mismatch.");
583        }
584    
585        // We need to consider the possibility that tags are missing or in the wrong order
586        // My guiding assumption here is: All tagged Datas are assumed to have the default value for
587        // all tags which are not explicitly defined
588    
589        const DataTagged* mptr=dynamic_cast<const DataTagged*>(mask2.m_data.get());
590        const DataTagged* optr=dynamic_cast<const DataTagged*>(other2.m_data.get());
591        DataTagged* tptr=dynamic_cast<DataTagged*>(m_data.get());
592    
593        // first, add any tags missing from other or mask
594        const DataTagged::DataMapType& olookup=optr->getTagLookup();
595            const DataTagged::DataMapType& mlookup=mptr->getTagLookup();
596        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
597        for (i=olookup.begin();i!=olookup.end();i++)
598        {
599               tptr->addTag(i->first);
600            }
601            for (i=mlookup.begin();i!=mlookup.end();i++) {
602               tptr->addTag(i->first);
603            }
604        // now we know that *this has all the required tags but they aren't guaranteed to be in
605        // the same order
606    
607        // There are two possibilities: 1. all objects have the same rank. 2. other is a scalar
608        if ((selfrank==otherrank) && (otherrank==maskrank))
609        {
610            for (i=olookup.begin();i!=olookup.end();i++)
611            {
612                // get the target offset
613                DataTypes::ValueType::size_type toff=tptr->getOffsetForTag(i->first);
614                    DataTypes::ValueType::size_type moff=mptr->getOffsetForTag(i->first);
615                DataTypes::ValueType::size_type ooff=optr->getOffsetForTag(i->first);
616                for (int i=0;i<getDataPointSize();++i)
617                {
618                    if (mvec[i+moff]>0)
619                    {
620                        self[i+toff]=ovec[i+ooff];
621                    }
622                }
623                }
624        }
625        else    // other is a scalar
626        {
627            for (i=mlookup.begin();i!=mlookup.end();i++)
628            {
629                // get the target offset
630                DataTypes::ValueType::size_type toff=tptr->getOffsetForTag(i->first);
631                    DataTypes::ValueType::size_type moff=mptr->getOffsetForTag(i->first);
632                DataTypes::ValueType::size_type ooff=optr->getOffsetForTag(i->first);
633                for (int i=0;i<getDataPointSize();++i)
634                {
635                    if (mvec[i+moff]>0)
636                    {
637                        self[i+toff]=ovec[ooff];
638                    }
639                }
640                }
641        }
642    
643        return;         // ugly
644      }
645      // mixed scalar and non-scalar operation
646      if ((selfrank>0) && (otherrank==0) && (mask2.getDataPointShape()==getDataPointShape()))
647      {
648            size_t num_points=self.size();
649        // OPENMP 3.0 allows unsigned loop vars.
650        #if defined(_OPENMP) && (_OPENMP < 200805)
651        long i;
652        #else
653        size_t i;
654        #endif      
655        #pragma omp parallel for private(i) schedule(static)
656        for (i=0;i<num_points;++i)
657        {
658            if (mvec[i]>0)
659            {
660                self[i]=ovec[0];
661            }
662        }
663        return;         // ugly!
664      }
665      // tagged data is already taken care of so we only need to worry about shapes
666      // special cases with scalars are already dealt with so all we need to worry about is shape
667      if ((getDataPointShape()!=other2.getDataPointShape()) || getDataPointShape()!=mask2.getDataPointShape())
668    {    {
669      throw DataException("Error - size mismatch in arguments to copyWithMask.");      ostringstream oss;
670        oss <<"Error - size mismatch in arguments to copyWithMask.";
671        oss << "\nself_shape=" << DataTypes::shapeToString(getDataPointShape());
672        oss << " other2_shape=" << DataTypes::shapeToString(other2.getDataPointShape());
673        oss << " mask2_shape=" << DataTypes::shapeToString(mask2.getDataPointShape());
674        throw DataException(oss.str());
675    }    }
676    size_t num_points=self.size();    size_t num_points=self.size();
677    

Legend:
Removed from v.2245  
changed lines
  Added in v.2246

  ViewVC Help
Powered by ViewVC 1.1.26