.. _program_listing_file_larflow_PrepFlowMatchData_FlowMatchHitMaker.cxx:

Program Listing for File FlowMatchHitMaker.cxx
==============================================

|exhale_lsh| :ref:`Return to documentation for file <file_larflow_PrepFlowMatchData_FlowMatchHitMaker.cxx>` (``larflow/PrepFlowMatchData/FlowMatchHitMaker.cxx``)

.. |exhale_lsh| unicode:: U+021B0 .. UPWARDS ARROW WITH TIP LEFTWARDS

.. code-block:: cpp

   #include "FlowMatchHitMaker.h"
   #define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
   #include <numpy/ndarrayobject.h>
   
   #include "larcv/core/PyUtil/PyUtils.h"
   #include "core/LArUtil/LArProperties.h"
   #include "core/LArUtil/Geometry.h"
   #include "larcv/core/Base/larcv_logger.h"
   
   #include <stdexcept>
   
   namespace larflow {
   namespace prep {
   
     int FlowMatchHitMaker::add_match_data( PyObject* pair_probs,
                                            PyObject* source_sparseimg, PyObject* target_sparseimg,
                                            PyObject* matchpairs,
                                            const int source_plane, const int target_plane,
                                            const larcv::ImageMeta& source_meta,
                                            const std::vector<larcv::Image2D>& img_v,
                                            const larcv::EventChStatus& ev_chstatus )
     {
   
       // enable numpy environment (if not already set)
       std::cout << "setting pyutils ... ";
       import_array1(0);    
       larcv::SetPyUtil();
       std::cout << " done" << std::endl;
   
       // cast numpy data to C-arrays
       const int dtype = NPY_FLOAT;
       PyArray_Descr *descr = PyArray_DescrFromType(dtype);
   
       // match scores
       npy_intp pair_dims[2];
       float **probs_carray;
       std::cout << "get pair prob dims" << std::endl;
       if ( PyArray_AsCArray( &pair_probs, (void**)&probs_carray, pair_dims, 2, descr )<0 ) {
         larcv::logger::get("PyFlowMatchData::make_larflow_hits").send(larcv::msg::kCRITICAL,__FUNCTION__,__LINE__, "cannot get carray for pair prob matrix");
       }
   
       // sparse source image
       npy_intp source_dims[2];
       float **source_carray;
       std::cout << "get source dims" << std::endl;    
       if ( PyArray_AsCArray( &source_sparseimg, (void**)&source_carray, source_dims, 2, descr )<0 ) {
         larcv::logger::get("PyFlowMatchData::make_larflow_hits").send(larcv::msg::kCRITICAL,__FUNCTION__,__LINE__, "cannot get carray for source sparse-image matrix");
       }
   
       // sparse target image
       npy_intp target_dims[2];
       float **target_carray;
       std::cout << "get target dims" << std::endl;        
       if ( PyArray_AsCArray( &target_sparseimg, (void**)&target_carray, target_dims, 2, descr )<0 ) {
         larcv::logger::get("PyFlowMatchData::make_larflow_hits").send(larcv::msg::kCRITICAL,__FUNCTION__,__LINE__, "cannot get carray for target sparse-image matrix");
       }
   
       // (source,target) indicies for each match
       npy_intp match_dims[2];
       long **matchpairs_carray;
       std::cout << "get index dims" << std::endl;            
       if ( PyArray_AsCArray( &matchpairs, (void**)&matchpairs_carray, match_dims, 2, PyArray_DescrFromType(NPY_LONG) )<0 ) {
         larcv::logger::get("PyFlowMatchData::make_larflow_hits").send(larcv::msg::kCRITICAL,__FUNCTION__,__LINE__, "cannot get carray for match pair matrix");
       }
   
       std::cout << "match matrix:  (" << match_dims[0]  << "," << match_dims[1]  << ")" << std::endl;    
       std::cout << "prob matrix:   (" << pair_dims[0]   << "," << pair_dims[1]   << ")" << std::endl;    
       std::cout << "source matrix: (" << source_dims[0] << "," << source_dims[1] << ")" << std::endl;
       std::cout << "target matrix: (" << target_dims[0] << "," << target_dims[1] << ")" << std::endl;
   
       // uboone-specific geo stuff
   
       int other_plane = -1;
       if ( source_plane==2 ) {
         other_plane = ( target_plane==0 ) ? 1 : 0;
       }
       else if ( source_plane==1 ) {
         other_plane = ( target_plane==0 ) ? 2 : 0;
       }
       else if ( source_plane==0 ) {
         other_plane = ( target_plane==1 ) ? 2 : 1;
       }
   
       auto const& other_plane_status = ev_chstatus.Status( other_plane );
   
       // choose threshold
       //float match_prob_threshold = 0.5;
       // if (other_plane==2) 
       //   match_prob_threshold = 0.1;
   
       bool require_3plane = false;
         
       // loop over each candidate match
       int n_indead = 0;
       int n_2plane = 0;
       int n_3plane = 0;
       int nrepeated = 0;
       for (int ipair=0; ipair<(int)pair_dims[1]; ipair++) {
   
         // score for this match
         float prob = probs_carray[0][ipair];
   
         // match threshold
         //if ( prob<_match_score_threshold ) continue;
   
         // get source and target index (in sparse image)
         int srcidx = matchpairs_carray[ipair][0];
         int taridx = matchpairs_carray[ipair][1];
   
         // get source and target column, row
         int srccol = (int)source_carray[srcidx][1];
         int srcrow = (int)source_carray[srcidx][0];      
         int tarcol = (int)target_carray[taridx][1];
   
         // get the other plane wire
         // convert row to tick
         double y, z;
         larutil::Geometry::GetME()->IntersectionPoint( srccol, tarcol, (UChar_t)source_plane, (UChar_t)target_plane, y, z );
         float tick = source_meta.pos_y( srcrow );
   
         if ( y<-117 || y>117 || z<0 || z>1036 ) continue;
   
         Double_t pos[3] = { 0, y, z };
         float other_wire = larutil::Geometry::GetME()->WireCoordinate( pos, other_plane );
         float other_adc  = 0;
         try {
           img_v[other_plane].pixel( srcrow, (int)other_wire, __FILE__, __LINE__ );
         }
         catch (std::exception& e ) {
           std::cout << "Error getting other plane pixel, (r,c)=(" << srcrow << "," << other_wire << "). "
                     << "for pos(y,z)=(" << y << ", " << z << "). "
                     << "from (src,tar) intersection=(" << srccol << "," << tarcol << ")"
                     << std::endl;
           std::cout << "error: " << e.what() << std::endl;
           continue;
         }
   
         std::vector<int> img_coords(3,0);
         img_coords[ source_plane ] = srccol;
         img_coords[ target_plane ] = tarcol;
         img_coords[ other_plane  ] = (int)other_wire;
   
         // std::cout << "pair[" << ipair << "] (t,y,z)=(" << tick << "," << y << "," << z << ") "
         //           << "imgcoords=[" << img_coords[0] << "," << img_coords[1] << "," << img_coords[2] << "]"
         //           << std::endl;
         //std::cout << "    otherplane=" << other_plane << std::endl;
   
         // determine if this is a valid point
         // (1) must have charge inall three planes OR
         // (2) charge in two and dead channel in one      
         bool indead = false;
         if ( other_adc<10.0 ) {
   
           // if event chstatus pointer is not null, we check the ch status
           int chstatus = other_plane_status.Status( other_wire );
           if ( chstatus==4 ) {
             if ( require_3plane ) {
               // good + below threshold so we ignore this match
               continue;
             }
             else {
               // accept anyway
               n_2plane++;
             }
           }
           else {
             // in dead, so accept
             indead = true;
             n_indead++;
           }
         }
         else {
           // all three planes have charge
           n_3plane++;
         }
   
         // now look for the wire-triple in the match data list
         std::vector<int> triple(4,0);
         triple[ source_plane ] = srccol;
         triple[ target_plane ] = tarcol;
         triple[ other_plane  ] = int(other_wire);
         triple[ 3 ] = (int)tick;
   
         //std::cout << "    triple=" << triple[0] << "," << triple[1] << "," << triple[2] << std::endl;
   
         auto it = _match_map.find( triple );
         if ( it==_match_map.end() ) {
           // if not in match map, we create a new entry
           match_t m;
           m.set_wire( source_plane, srccol );
           m.set_wire( target_plane, tarcol );
           m.set_wire( other_plane,  int(other_wire) );
           m.tyz = { (float)tick, (float)y, (float)z };
           _matches_v.emplace_back( std::move(m) );
           _match_map[ triple ] = int(_matches_v.size())-1;
           it = _match_map.find( triple );
         }
         else {
           nrepeated++;
         }
   
         auto& m = _matches_v.at( it->second );
         // set the score given here
         //if ( indead ) prob *= 0.5;
         m.set_score( source_plane, target_plane, prob );
   
       }//end of score loop
   
   
       PyArray_Free( pair_probs, (void*)probs_carray );
       PyArray_Free( source_sparseimg, (void*)source_carray );
       PyArray_Free( target_sparseimg, (void*)target_carray );
       PyArray_Free( matchpairs, (void*)matchpairs_carray );    
       
       std::cout << "-------------------------------------------------" << std::endl;
       std::cout << "process match data. " << std::endl;
       std::cout << "  number of triples with match information: " << _match_map.size() << std::endl;
       std::cout << "  number of: 3-plane=" << n_3plane << " 2-plane=" << n_2plane << " 2-plane+dead=" << n_indead << std::endl;
       std::cout << "  num repeated: " << nrepeated << std::endl;
       std::cout << "-------------------------------------------------" << std::endl;
     }
   
     void FlowMatchHitMaker::make_hits( const larcv::EventChStatus& ev_chstatus,
                                        std::vector<larlite::larflow3dhit>& hit_v ) const {
   
       const float cm_per_tick = larutil::LArProperties::GetME()->DriftVelocity()*0.5;
       const int ncolumns = 19;
       
       int idx = 0;
       unsigned long maxsize = hit_v.size() + _matches_v.size()+10;
       hit_v.reserve(maxsize);
       for ( auto const& m : _matches_v ) {
   
         // find the highest match score
         std::vector<float> scores = m.get_scores();
         float maxscore = 0;
         int maxdir = -1;
         for ( int i=0; i<scores.size(); i++ ) {
           float s = scores[i];
           if ( s>maxscore ) {
             maxscore = s;
             maxdir = i;
           }
         }
         if ( maxscore<_match_score_threshold )
           continue;
         
         larlite::larflow3dhit hit;
         hit.tick = m.tyz[0];
         hit.srcwire = m.Y;
         hit.targetwire.resize(3);
         hit.targetwire[0] = m.U;
         hit.targetwire[1] = m.V;
         hit.targetwire[2] = m.Y;
         hit.idxhit = idx;
         if(m.istruth==0) hit.truthflag = larlite::larflow3dhit::TruthFlag_t::kNoTruthMatch;
         else{ hit.truthflag = larlite::larflow3dhit::TruthFlag_t::kOnTrack;}
         float x = (m.tyz[0]-3200.0)*cm_per_tick;
         hit.resize(ncolumns,0);
         hit[0] = x;
         hit[1] = m.tyz[1];
         hit[2] = m.tyz[2];
         
         // store scores
         for ( int i=0; i<scores.size(); i++ ) {
           float s = scores[i];
           hit[3+i] = s;
         }
         hit.flowdir = (larlite::larflow3dhit::FlowDirection_t)maxdir;
         hit[9] = maxscore;
         hit.track_score = maxscore;
   
         if ( has_ssnet_scores && m.ssnet_scores.size()==3) {
           for (int c=0; c<3; c++)
             hit[10 + c] = m.ssnet_scores[c];
         }
         
         if ( has_kplabel_scores ) {
           for (int c=0; c<3; c++) 
             hit[13+c] = m.keypoint_scores[c];
         }
   
         if ( has_paf ) {
           for (int i=0; i<3; i++)
             hit[16+i] = m.paf[i];
         }
         
         hit_v.emplace_back( std::move(hit) );
       }
       std::cout << "[FlowMatchHitMaker::make_hits] saved " << hit_v.size() << " hits "
                 << " from "  << _matches_v.size() << " matches" << std::endl;
       
     };
   
     int FlowMatchHitMaker::add_triplet_match_data( PyObject* triple_probs,
                                                    PyObject* triplet_indices,
                                                    PyObject* imgu_sparseimg,
                                                    PyObject* imgv_sparseimg,
                                                    PyObject* imgy_sparseimg,
                                                    const std::vector< std::vector<float> >& pos_vv,
                                                    const std::vector<larcv::Image2D>& adc_v ) {
   
       // enable numpy environment (if not already set)
       std::cout << "setting pyutils ... ";
       import_array1(0);    
       larcv::SetPyUtil();
       std::cout << " done" << std::endl;
   
       // // match scores
       int pair_ndims = PyArray_NDIM( (PyArrayObject*)triplet_indices );
       npy_intp* pair_dims = PyArray_DIMS( (PyArrayObject*)triplet_indices );
       std::cout << "[FlowMatchHitMaker] pair prob dims=(" << pair_dims[0] << ")" << std::endl;
   
       // // triplet indicies for each match
       // npy_intp index_dims[2];
       // std::cout << "[FlowMatchHitMaker] index array dims=(" << index_dims[0] << "," << index_dims[1] << ")" << std::endl;
       
       // // sparse images
       // npy_intp imgu_dims[2];
       // std::cout << "[FlowMatchHitMaker] img[u] prob dims=(" << imgu_dims[0] << "," << imgu_dims[1] << ")" << std::endl;
   
       // npy_intp imgv_dims[2];
       // std::cout << "[FlowMatchHitMaker] img[v] prob dims=(" << imgv_dims[0] << "," << imgv_dims[1] << ")" << std::endl;
   
       // npy_intp imgy_dims[2];
       // std::cout << "[FlowMatchHitMaker] img[y] prob dims=(" << imgy_dims[0] << "," << imgy_dims[1] << ")" << std::endl;
   
       bool precalc_pos = ( pos_vv.size()>0 ) ? true : false;
       int nrepeated = 0;
       int nbelowminprob = 0.;
       
       // loop over each candidate triplet
       for (int ipair=0; ipair<(int)pair_dims[0]; ipair++) {
   
         // score for this match
         float prob = *(float*)PyArray_GETPTR1( (PyArrayObject*)triple_probs, ipair );
   
         long index[4] = { *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 0 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 1 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 2 ),
               *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 3 ) };
   
         
         std::vector<int> triple(5,0); // (col,col,col,tick,istruth)
         triple[ 0 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgu_sparseimg, index[0], 1 );
         triple[ 1 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgv_sparseimg, index[1], 1 );
         triple[ 2 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 1 );
         int row = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 0 );
         triple[ 3 ] = (int)adc_v[2].meta().pos_y( row );
         triple[ 4 ] = (int)index[3]; //truth 
         
         // match threshold
         if ( prob<_match_score_threshold ) {
           nbelowminprob++;
           //continue;
         }
   
         std::vector<float> pos(3,0);
         
         if ( precalc_pos ) {
           pos = pos_vv[ipair];
         }
         else {
           
           // get the other plane wire
           // convert row to tick
           double y, z;
           larutil::Geometry::GetME()->IntersectionPoint( triple[2], triple[0], (UChar_t)2, (UChar_t)0, y, z );
   
           pos[0] = (triple[3]-3200.0)*0.5*larutil::LArProperties::GetME()->DriftVelocity();
           pos[1] = y;
           pos[2] = z;
           
         }
   
         auto it = _match_map.find( triple );
         if ( it==_match_map.end() ) {
           // if not in match map, we create a new entry
           match_t m;
       m.set_istruth(triple[4]);
           for ( size_t p=0; p<3; p++ )
             m.set_wire( p, triple[p] );
           m.tyz = { (float)triple[3], (float)pos[1], (float)pos[2] };
           _matches_v.emplace_back( std::move(m) );
           _match_map[ triple ] = int(_matches_v.size())-1;
           it = _match_map.find( triple );
         }
         else {
           std::cout << "repeated triple: (" << triple[0] << "," << triple[1] << "," << triple[2] << "," << triple[3] << "." << triple[4] <<")" << std::endl;
           nrepeated++;
         }
   
         auto& m = _matches_v.at( it->second );
         // set the score given here
         //if ( indead ) prob *= 0.5;
         m.set_score( 2, 0, prob );
   
       }//end of score loop
       
       std::cout << "-------------------------------------------------" << std::endl;
       std::cout << "process match data. " << std::endl;
       std::cout << "  number of triples with match information: " << _match_map.size() << std::endl;
       std::cout << "  num repeated: " << nrepeated << std::endl;
       std::cout << "  num below min prob: " << nbelowminprob << std::endl;
       std::cout << "-------------------------------------------------" << std::endl;
     }
   
   
     int FlowMatchHitMaker::add_triplet_ssnet_scores( PyObject* triplet_indices,
                                                      PyObject* imgu_sparseimg,
                                                      PyObject* imgv_sparseimg,
                                                      PyObject* imgy_sparseimg,
                                                      const larcv::ImageMeta& meta,
                                                      PyObject* ssnet_scores ) {
   
       has_ssnet_scores = true;
       
       // match triplets
       int pair_ndims      = PyArray_NDIM( (PyArrayObject*)triplet_indices );
       npy_intp* pair_dims = PyArray_DIMS( (PyArrayObject*)triplet_indices );
   
       int ss_pair_ndims      = PyArray_NDIM( (PyArrayObject*)ssnet_scores );
       npy_intp* ss_pair_dims = PyArray_DIMS( (PyArrayObject*)ssnet_scores );
       std::cout << "[FlowMatchHitMaker::add_triplet_ssnet_scores] score dims="
                 << "(" << ss_pair_dims[0] << "," << ss_pair_dims[1] << ")"
                 << std::endl;            
   
       for (int ipair=0; ipair<(int)ss_pair_dims[1]; ipair++ ) {
         // get the triplet indices
         
         long index[4] = { *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 0 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 1 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 2 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 3 ) };
         
         std::vector<int> triple(5,0); // (col,col,col,tick,istruth)
         triple[ 0 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgu_sparseimg, index[0], 1 );
         triple[ 1 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgv_sparseimg, index[1], 1 );
         triple[ 2 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 1 );
         int row = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 0 );
         triple[ 3 ] = (int)meta.pos_y( row );
         triple[ 4 ] = (int)index[3];
       
         // look for the triplet in the map
         auto it = _match_map.find( triple );
         if ( it!=_match_map.end() ) {
           int matchidx = it->second;
           auto& match = _matches_v[matchidx];
   
           // store the ssnet scores
           match.ssnet_scores.resize(3,0.0);
           //std::cout << "save ssnet scores: ";
           for (int p=0; p<3; p++) {
             match.ssnet_scores[p] = (float) *((float*)PyArray_GETPTR2( (PyArrayObject*)ssnet_scores, p, ipair ));
             //std::cout << match.ssnet_scores[p] << " ";
           }
           //std::cout << std::endl;
         }
         else {
           std::cout << "could not find match for triplet: ("
                     << triple[0] << ","
                     << triple[1] << ","
                     << triple[2] << ","
                     << triple[3] << ","
                     << triple[4] << ")"
                     << std::endl;
           throw std::runtime_error("bad ssnet data");
         }
       }
       
       return 0;
     }
   
     int FlowMatchHitMaker::add_triplet_keypoint_scores( PyObject* triplet_indices,
                                                         PyObject* imgu_sparseimg,
                                                         PyObject* imgv_sparseimg,
                                                         PyObject* imgy_sparseimg,
                                                         const larcv::ImageMeta& meta,                                                      
                                                         PyObject* kplabel_scores ) {
   
       has_kplabel_scores = true;
       
       // match triplets
       int pair_ndims      = PyArray_NDIM( (PyArrayObject*)triplet_indices );
       npy_intp* pair_dims = PyArray_DIMS( (PyArrayObject*)triplet_indices );
   
       int kp_pair_ndims      = PyArray_NDIM( (PyArrayObject*)kplabel_scores );
       npy_intp* kp_pair_dims = PyArray_DIMS( (PyArrayObject*)kplabel_scores );
       std::cout << "[FlowMatchHitMaker::add_triplet_keypoint_scores] score dims=(" << kp_pair_dims[0] << ")" << std::endl;    
   
       for (int ipair=0; ipair<(int)kp_pair_dims[0]; ipair++ ) {
         // get the triplet indices
         
         long index[4] = { *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 0 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 1 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 2 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 3 ) };
         
         std::vector<int> triple(5,0); // (col,col,col,tick,istruth)
         triple[ 0 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgu_sparseimg, index[0], 1 );
         triple[ 1 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgv_sparseimg, index[1], 1 );
         triple[ 2 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 1 );
         int row = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 0 );
         triple[ 3 ] = (int)meta.pos_y( row );
         triple[ 4 ] = (int)index[3];
       
         // look for the triplet in the map
         auto it = _match_map.find( triple );
         if ( it!=_match_map.end() ) {
           int matchidx = it->second;
           auto& match = _matches_v[matchidx];
   
           // store the kplabel scores]
           if ( pair_ndims>1 ) {
             match.keypoint_scores.resize(pair_dims[1],0);
             for (int iclass=0; iclass<(int)pair_dims[1]; iclass++)
               match.keypoint_scores[iclass] = (float)*((float*)PyArray_GETPTR2( (PyArrayObject*)kplabel_scores, ipair, iclass ));
           }
           else {
             match.keypoint_scores.resize(1,0);
             match.keypoint_scores[0] = (float)*((float*)PyArray_GETPTR1( (PyArrayObject*)kplabel_scores, ipair ) );
           }
         }
         else {
           std::cout << "could not find match for triplet: ("
                     << triple[0] << ","
                     << triple[1] << ","
                     << triple[2] << ","
                     << triple[3] << ","
                     << triple[4] << ")"
                     << std::endl;
         }
       }
       
       return 0;
     }
   
     int FlowMatchHitMaker::add_triplet_affinity_field( PyObject* triplet_indices,
                                                        PyObject* imgu_sparseimg,
                                                        PyObject* imgv_sparseimg,
                                                        PyObject* imgy_sparseimg,
                                                        const larcv::ImageMeta& meta,                                                      
                                                        PyObject* paf_pred ) {
   
       has_paf = true;
       
       // match triplets
       int pair_ndims      = PyArray_NDIM( (PyArrayObject*)triplet_indices );
       npy_intp* pair_dims = PyArray_DIMS( (PyArrayObject*)triplet_indices );
   
       int paf_pair_ndims      = PyArray_NDIM( (PyArrayObject*)paf_pred );
       npy_intp* paf_pair_dims = PyArray_DIMS( (PyArrayObject*)paf_pred );
       std::cout << "[FlowMatchHitMaker::add_triplet_affinity_field] score dims=(" << paf_pair_dims[0] << ")" << std::endl;    
   
       for (int ipair=0; ipair<(int)paf_pair_dims[0]; ipair++ ) {
         // get the triplet indices
         
         long index[4] = { *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 0 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 1 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 2 ),
                           *(long*)PyArray_GETPTR2( (PyArrayObject*)triplet_indices, ipair, 3 ) };
         
         std::vector<int> triple(5,0); // (col,col,col,tick,istruth)
         triple[ 0 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgu_sparseimg, index[0], 1 );
         triple[ 1 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgv_sparseimg, index[1], 1 );
         triple[ 2 ] = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 1 );
         int row = (int)*(float*)PyArray_GETPTR2( (PyArrayObject*)imgy_sparseimg, index[2], 0 );
         triple[ 3 ] = (int)meta.pos_y( row );
         triple[ 4 ] = (int)index[3];
       
         // look for the triplet in the map
         auto it = _match_map.find( triple );
         if ( it!=_match_map.end() ) {
           int matchidx = it->second;
           auto& match = _matches_v[matchidx];
   
           // store the kplabel scores]
           match.paf.resize(3,0);
           for (int i=0;i<3; i++ ) {
             match.paf[i] = (float)*((float*)PyArray_GETPTR2( (PyArrayObject*)paf_pred, ipair, i ));
           }
         }
         else {
           std::cout << "could not find match for triplet: ("
                     << triple[0] << ","
                     << triple[1] << ","
                     << triple[2] << ","
                     << triple[3] << ","
                     << triple[4] << ")"
                     << std::endl;
         }
       }
       
       return 0;
     }
     
   }  
   }