Changeset 2567 for branches


Ignore:
Timestamp:
06/14/12 10:01:40 (12 years ago)
Author:
Takeshi Nakazato
Message:

New Development: No

JIRA Issue: No

Ready for Test: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Avoid division by zero.


File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/hpc33/src/STMath.cpp

    r2566 r2567  
    47814781  const uInt *p = rows.data() ;
    47824782  vector<int> ids( 2 ) ;
     4783  Block<uInt> flagchan( spsize ) ;
     4784  uInt nflag = 0 ;
    47834785  for ( int irow = 0 ; irow < rows.nelements() ; irow++ ) {
    47844786    double reftime = timeCol.asdouble(*p) ;
    … …  
    47954797      // using gain array
    47964798      for ( unsigned int j = 0 ; j < tcal.size() ; j++ ) {
    4797         spec[j] = ( ( spec[j] - spoff[j] ) / spoff[j] )
    4798           * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     4799        if ( spoff[j] == 0.0 || (sphot[j]-spsky[j]) == 0.0 ) {
     4800          spec[j] = 0.0 ;
     4801          flagchan[nflag++] = j ;
     4802        }
     4803        else {
     4804          spec[j] = ( ( spec[j] - spoff[j] ) / spoff[j] )
     4805            * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     4806        }
    47994807      }
    48004808    }
    … …  
    48024810      // Chopper-Wheel calibration (Ulich & Haas 1976)
    48034811      for ( unsigned int j = 0 ; j < tcal.size() ; j++ ) {
    4804         spec[j] = ( spec[j] - spoff[j] ) / ( sphot[j] - spsky[j] ) * tcal[j] ;
     4812        if ( (sphot[j]-spsky[j]) == 0.0 ) {
     4813          spec[j] = 0.0 ;
     4814          flagchan[nflag++] = j ;
     4815        }
     4816        else {
     4817          spec[j] = ( spec[j] - spoff[j] ) / ( sphot[j] - spsky[j] ) * tcal[j] ;
     4818        }
    48054819      }
    48064820    }
    48074821    out->specCol_.put( *p, spec ) ;
    48084822    out->tsysCol_.put( *p, tsys ) ;
     4823    if ( nflag > 0 ) {
     4824      Vector<uChar> fl = out->flagsCol_( *p ) ;
     4825      for ( unsigned int j = 0 ; j < nflag ; j++ ) {
     4826        fl[flagchan[j]] = (uChar)True ;
     4827      }
     4828      out->flagsCol_.put( *p, fl ) ;
     4829    }
     4830    nflag = 0 ;
    48094831    p++ ;
    48104832  }
    … …  
    48314853  const uInt *p = rows.data() ;
    48324854  vector<int> ids( 2 ) ;
     4855  Block<uInt> flagchan( spsize ) ;
     4856  uInt nflag = 0 ;
    48334857  for ( int irow = 0 ; irow < rows.nelements() ; irow++ ) {
    48344858    double reftime = timeCol.asdouble(*p) ;
    … …  
    48454869    unsigned int tsyssize = tsys.nelements() ;
    48464870    for ( unsigned int j = 0 ; j < spsize ; j++ ) {
    4847       spec[j] = ( spec[j] - spoff[j] ) / spoff[j] ;
     4871      if ( spoff[j] == 0.0 ) {
     4872        spec[j] = 0.0 ;
     4873        flagchan[nflag++] = j ;
     4874      }
     4875      else {
     4876        spec[j] = ( spec[j] - spoff[j] ) / spoff[j] ;
     4877      }
    48484878      if ( tsyssize == spsize )
    48494879        spec[j] *= tsys[j] ;
    … …  
    48524882    }
    48534883    out->specCol_.put( *p, spec ) ;
     4884    if ( nflag > 0 ) {
     4885      Vector<uChar> fl = out->flagsCol_( *p ) ;
     4886      for ( unsigned int j = 0 ; j < nflag ; j++ ) {
     4887        fl[flagchan[j]] = (uChar)True ;
     4888      }
     4889      out->flagsCol_.put( *p, fl ) ;
     4890    }
     4891    nflag = 0 ;
    48544892    p++ ;
    48554893  }
    … …  
    48904928  const uInt *p = rows.data() ;
    48914929  vector<int> ids( 2 ) ;
     4930  Block<uInt> flagchan( spsize ) ;
     4931  uInt nflag = 0 ;
    48924932  for ( int irow = 0 ; irow < rows.nelements() ; irow++ ) {
    48934933    double reftime = timeCol.asdouble(*p) ;
    … …  
    49084948    Vector<Float> spref = on[1]->specCol_( *p ) ;
    49094949    for ( unsigned int j = 0 ; j < spsize ; j++ ) {
    4910       spec[j] = tcalS[j] * spsig[j] / ( sphotS[j] - spskyS[j] )
    4911         - tcalR[j] * spref[j] / ( sphotR[j] - spskyR[j] ) ;
     4950      if ( (sphotS[j]-spskyS[j]) == 0.0 || (sphotR[j]-spskyR[j]) == 0.0 ) {
     4951        spec[j] = 0.0 ;
     4952        flagchan[nflag++] = j ;
     4953      }
     4954      else {
     4955        spec[j] = tcalS[j] * spsig[j] / ( sphotS[j] - spskyS[j] )
     4956          - tcalR[j] * spref[j] / ( sphotR[j] - spskyR[j] ) ;
     4957      }
    49124958    }
    49134959    sig->specCol_.put( *p, spec ) ;
    … …  
    49164962    ref->specCol_.put( *p, spec ) ;
    49174963    ref->tsysCol_.put( *p, tsysR ) ;   
     4964    if ( nflag > 0 ) {
     4965      Vector<uChar> flsig = sig->flagsCol_( *p ) ;
     4966      Vector<uChar> flref = ref->flagsCol_( *p ) ;
     4967      for ( unsigned int j = 0 ; j < nflag ; j++ ) {
     4968        flsig[flagchan[j]] = (uChar)True ;
     4969        flref[flagchan[j]] = (uChar)True ;
     4970      }
     4971      sig->flagsCol_.put( *p, flsig ) ;
     4972      ref->flagsCol_.put( *p, flref ) ;
     4973    }
     4974    nflag = 0 ;
    49184975    p++ ;
    49194976  }
    … …  
    49475004  const uInt *p = rows.data() ;
    49485005  vector<int> ids( 2 ) ;
     5006  Block<uInt> flagchan( spsize ) ;
     5007  uInt nflag = 0 ;
    49495008  for ( int irow = 0 ; irow < rows.nelements() ; irow++ ) {
    49505009    double reftime = timeCol.asdouble(*p) ;
    … …  
    49595018    // using gain array
    49605019    for ( unsigned int j = 0 ; j < spsize ; j++ ) {
    4961       spec[j] = ( ( spsig[j] - spref[j] ) / spref[j] )
    4962         * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     5020      if ( spref[j] == 0.0 || (sphot[j]-spsky[j]) == 0.0 ) {
     5021        spec[j] = 0.0 ;
     5022        flagchan[nflag++] = j ;
     5023      }
     5024      else {
     5025        spec[j] = ( ( spsig[j] - spref[j] ) / spref[j] )
     5026          * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     5027      }
    49635028    }
    49645029    sig->specCol_.put( *p, spec ) ;
    49655030    sig->tsysCol_.put( *p, tsys ) ;
     5031    if ( nflag > 0 ) {
     5032      Vector<uChar> fl = sig->flagsCol_( *p ) ;
     5033      for ( unsigned int j = 0 ; j < nflag ; j++ ) {
     5034        fl[flagchan[j]] = (uChar)True ;
     5035      }
     5036      sig->flagsCol_.put( *p, fl ) ;
     5037    }
     5038    nflag = 0 ;
    49665039
    49675040    reftime = timeCol2.asdouble(*p) ;
    … …  
    49735046    // using gain array
    49745047    for ( unsigned int j = 0 ; j < spsize ; j++ ) {
    4975       spec[j] = ( ( spref[j] - spsig[j] ) / spsig[j] )
    4976         * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     5048      if ( spsig[j] == 0.0 || (sphot[j]-spsky[j]) == 0.0 ) {
     5049        spec[j] = 0.0 ;
     5050        flagchan[nflag++] = j ;
     5051      }
     5052      else {
     5053        spec[j] = ( ( spref[j] - spsig[j] ) / spsig[j] )
     5054          * ( spsky[j] / ( sphot[j] - spsky[j] ) ) * tcal[j] ;
     5055      }
    49775056    }
    49785057    ref->specCol_.put( *p, spec ) ;
    49795058    ref->tsysCol_.put( *p, tsys ) ;   
     5059    if ( nflag > 0 ) {
     5060      Vector<uChar> fl = ref->flagsCol_( *p ) ;
     5061      for ( unsigned int j = 0 ; j < nflag ; j++ ) {
     5062        fl[flagchan[j]] = (uChar)True ;
     5063      }
     5064      ref->flagsCol_.put( *p, fl ) ;
     5065    }
     5066    nflag = 0 ;
    49805067    p++ ;
    49815068  }
Note: See TracChangeset for help on using the changeset viewer.