dcsimg
CodeGuru Home VC++ / MFC / C++ .NET / C# Visual Basic VB Forums Developer.com
Results 1 to 5 of 5

Thread: Azimuth search

  1. #1
    Join Date
    Aug 2012
    Posts
    5

    Azimuth search

    Dear friends.

    Good morning. I need little help.
    I am working to make full azimuth search. Can you help me in changing the code for this. I have attached file. Please help me.

    Looking forward to have help from you ..

    Code:
    //////////////////////////////////////////////////////////////
    
    #include"main.h"
    
    #define mymin(a,b) (a < b ? a : b) 
    #define mymax(a,b) (a > b ? a : b) 
    #define mysqr(a)  (a*a)
    
    #define  PI 3.1415927
    
    using namespace std;
    
    int main (int argc, char** argv) {
    
    ///MPI begin
      MPI_Init(&argc,&argv);
      int mpi_rank=0;
      MPI_Comm_rank(MPI_COMM_WORLD,&mpi_rank);
      int mpi_size;
      MPI_Comm_size(MPI_COMM_WORLD,&mpi_size);
    
      if(argc==1){
        PipeToMore(doc);
        return 0;
      }
      
    ///declarationen of variables///////////////
    
       float Gx, Gy;    
       float Sx, Sy,s;
       double starttime = 0.0;
       float velmin, velmax, *vel, dv, v;
       double endtime = 0.0;
       float dt, dx, apert;
       float x, xm,tt, t0;
       int ns, nv, s_count, rcoher, time, flag;
       int idown, iup;
       string str;
       string Inputfile,  predictionfile;
       string Outputfile;
       float Value, Val_up, Val_down, Value_temp, Val_temp;
       float *sembl_up, *sembl_down;
       string Output_name[mpi_size];
       float apery, aperx;
    
    ///Getting parameter /////////////////////////////////////////////////////////////////////
      Parser* parser = new Parser(argc, argv);
      int *intdummy;
      float *fdummy;
    
    
       Inputfile = parser->parse("input");
        if (!Inputfile.length()) {
             if (mpi_rank==0)
        cerr << "ERROR: No In-filename found." << endl;
          MPI_Abort(MPI_COMM_WORLD,1);
        }
    
       if (!parser->parse("coherR",intdummy)) {
        if (mpi_rank==0) cerr << "ERROR: No coherence band width  found. Use coherR=" << endl;
        MPI_Abort(MPI_COMM_WORLD,1);
       }
       rcoher=*intdummy;
       
      if (!parser->parse("mode",intdummy)) {
        if (mpi_rank==0) cerr << "ERROR: No mode found. Use mode= 1 for 3D od 0 for 2D" << endl;
        MPI_Abort(MPI_COMM_WORLD,1);
       }
       int moda=*intdummy;
    
    
     if (!parser->parse("aperx",fdummy)) {
          if (mpi_rank==0)  cerr << "ERROR: No X-aperture found. Use aperx =" << endl;
        MPI_Abort(MPI_COMM_WORLD,1);
       }
       aperx=apert=*fdummy;
    
    if (moda ==1){
     if (!parser->parse("apery",fdummy)) {
          if (mpi_rank==0)  cerr << "ERROR: No Y-aperture found. Use apery =" << endl;
        MPI_Abort(MPI_COMM_WORLD,1);
       }
    
       apery=*fdummy;
    }
      
    if (!parser->parse("cmpdx",fdummy)) {
          if (mpi_rank==0)  cerr << "ERROR: No cmp distance found. Use cmpdx =" << endl;
        MPI_Abort(MPI_COMM_WORLD,1);
       }
      dx=*fdummy;
    
    
       Outputfile = parser->parse("output");
        if (!Outputfile.length()) {
         if (mpi_rank==0) 
      cerr << "ERROR: No output-filename found." << endl;
         MPI_Abort(MPI_COMM_WORLD,1);
       }
      
    float *vtarget;
       if(!parser->parse("vmin",vtarget)) 
      {
         if (mpi_rank==0) 
       cerr << "ERROR: No velocity found." << endl;
         MPI_Abort(MPI_COMM_WORLD,1);
       }
       velmin=*vtarget;
     if(!parser->parse("vmax",vtarget)) 
      {
         if (mpi_rank==0) 
           cerr << "ERROR: No velocity found." << endl;
         MPI_Abort(MPI_COMM_WORLD,1);
       }
       velmax=*vtarget;
    
      float vel_temp;
      if( velmax < velmin) 
      {
          vel_temp=velmax; velmax=velmin; velmin=vel_temp;
       }
    
        if(!parser->parse("dv",vtarget)) 
      {
         if (mpi_rank==0) 
    cerr << "ERROR: No velocity found." << endl;
         MPI_Abort(MPI_COMM_WORLD,1);
       }
      dv=*vtarget;
      
       delete parser;
       parser=0;
      
    
    /// SU data read in...
       SUdata *Input_data, *Output;
       Input_data = new SUdata(Inputfile.c_str());
       Input_data->readFile();
       SUdata *Output_data[mpi_size];
    
    if ( moda==0 ) apert=aperx;
    
    ///allocate all arrays and variable's assignment
       vector<Trace>::iterator Input_it, Input_temp, it_left, it_right;
       dt = float(Input_data->dt);
       ns = int(Input_data->nt);
       short fac1=short(Input_data->scalco);
       nv= static_cast<int>((velmax -velmin)/dv);
       vel = new float[nv+1];
       float dtnew=dt/1000000.;
       int aper = static_cast<int>(apert/dx);
       int width =  rcoher*2 ;
       sembl_up = new float[width+1];
       sembl_down = new float[width+1];
       float fac, ym, y, t_aper;
    
       if ( ( fac1 )  || ( fac1 != 0 ) )  fac=1./(1.*abs(fac1));
       else fac=1.;
    ///create string array for output name
      for (int i=mpi_rank; i<mpi_size;i++)
       {
         Output_name[i]+=Outputfile.c_str();
         str=ItoA(i);
         Output_name[i]+=str;
        }
    
    
     if (mpi_rank==0) {
         cout  << endl
              << "INFO: Velocity model building" << endl
           << "INFO: Copyright (C)  Stefan Duemmong & Sergius Dell - 2009" << endl
               << "INFO: sergius.dell@zmaw.de" << endl
              << endl;
        cout   << endl<<"Starting program on "<<mpi_size<<" nodes."<<endl<<endl;
        cout   << endl
           << "INFO: Using: " <<  endl
               << "INFO: kind of processing:     time-migration velocity analysis"  <<endl
               << "INFO: infile:          " << Inputfile << endl
              << "INFO: outfile:       " << Outputfile << endl
               << "INFO: cmp distance:           "  <<dx<<endl
               << "INFO: scaling:           "  <<fac<<endl
               << "INFO: radius of coherence band:           "  <<rcoher<<endl
               << "INFO: X aperture:           "  <<2*aperx<<endl;
    if (moda==1) cout<< "INFO: Y aperture:           "  <<2*apery<<endl
               << "INFO: mimimal velocity:           "  <<velmin<<endl
               << "INFO: maximal velocity:           "  <<velmax<<endl
               << "INFO: search step:           "  <<dv<<endl
               << endl;
      }
    ///allocate output data
    for (int i=mpi_rank; i<mpi_size;i++)
       {
        Output_data[i]=new SUdata(Output_name[i].c_str());
       }
    
    
     for (int i=0;i<=nv;i++) { 
       vel[i]=velmin + i*int(dv); 
     }
    
    ///Time start
       starttime=MPI_Wtime();
     
    
    if (moda==0 ) {
    ///Loop ovel all velocities
    
     for (int j=0; j<=nv;j++) {
         v=vel[j];
       if (mpi_rank==0)
               cout << "INFO: 2D semblance analysis:   velocity " << v<<" is done"<<endl;
    /// Loop over all traces: //////////////////////////////////////////////////
       for ( int i=mpi_rank; i<Input_data->_traces.size(); i+=mpi_size) {
    
         Input_it = Input_data->_traces.begin() + i;
    
          /// Output preparing: ///////////////////////////////////////////////////////////
               Gx = Input_it->gx;
           Sx = Input_it->sx;
           xm=0.5*(Gx+Sx)*fac;
             Output_data[mpi_rank]->_traces.push_back( Trace(ns));
           Output_data[mpi_rank]->_traces.back().gx = int(Gx);
               Output_data[mpi_rank]->_traces.back().sx = int(Sx);
           Output_data[mpi_rank]->_traces.back().dt = int(dt);
               Output_data[mpi_rank]->_traces.back().offset = int(v);
               Output_data[mpi_rank]->_traces.back().cdp = int(Input_it->cdp);    
               Output_data[mpi_rank]->_traces.back().scalco = short(fac1);
               ( i < aper )? it_left=Input_data->_traces.begin() +1 :it_left=Input_it - aper;
               ( (Input_data->_traces.size() - i) > aper )? it_right=Input_it + aper:it_right=Input_data->_traces.end()-1;
    
    /// Loop over all samples: /////////////////////////////////////////////////
           for (int ii=0; ii<ns; ii++)
               {
               t0=ii*dtnew;
               Value=0.;
                   s_count=0;
            for (int is=0; is < width;is++){
                  sembl_up[is]=0.;sembl_down[is]=0.;
                }
            for ( Input_temp=it_left;Input_temp != it_right;Input_temp++) {
                            x=0.5*(Input_temp->gx  + Input_temp->sx)*fac ;                                          /// xm difference
                tt=sqrt(t0*t0+4.*(xm-x)*(xm-x)/(v*v));                                              /// calculate tt hyperbola over all times 
                         idown=static_cast<int>(tt/dtnew);                                                   ///calculate the current index in arrays
                if (Input_temp > Input_data->_traces.end())  break;
                            ( (idown -rcoher ) > 0)? idown=idown-rcoher: idown=idown; 
                            ( (idown + rcoher) > ns)? idown=ns - rcoher: idown=idown;
                s_count++;                                                                /// trace counter
                for (int is=0; is<width;is++){                                                     /// semblance within the coherence band
                                Value_temp=float(Input_temp->data[idown+is]); 
                              sembl_up[is]+=Value_temp;
                              sembl_down[is]+=Value_temp*Value_temp;
                }
                            s=1.*s_count; 
                         }
             Val_up= Val_down=0.0;
                   for (int is=0; is<width;is++){
                 Val_up+=sembl_up[is]*sembl_up[is];
                  Val_down+=sembl_down[is];
             }
              Value=Val_up/Val_down/s; 
                  Output_data[mpi_rank]->_traces.back().data[ii]=Value; 
               }
       }
     }
    
      }
    if (moda==1 ) {
    ///Loop ovel all velocities
    
     for (int j=0; j<=nv;j++) {
         v=vel[j];
       if (mpi_rank==0)
               cout << "INFO: 3D semblance analysis:   velocity " << v<<" is done"<<endl;
    /// Loop over all traces: //////////////////////////////////////////////////
       for ( int i=mpi_rank; i<Input_data->_traces.size(); i+=mpi_size) {
         Input_it = Input_data->_traces.begin() + i;
    
          /// Output preparing: ///////////////////////////////////////////////////////////
               Gx = Input_it->gx;
           Sx = Input_it->sx;
               Gy = Input_it->gy;
           Sy = Input_it->sy;
           xm=0.5*(Gx+Sx)*fac;
               ym=0.5*(Gy+Sy)*fac;
           Output_data[mpi_rank]->_traces.push_back( Trace(ns));
           Output_data[mpi_rank]->_traces.back().gx = int(Gx);
               Output_data[mpi_rank]->_traces.back().sx = int(Sx);
           Output_data[mpi_rank]->_traces.back().gy = int(Gy);
               Output_data[mpi_rank]->_traces.back().sy = int(Sy);
           Output_data[mpi_rank]->_traces.back().dt = int(dt);
               Output_data[mpi_rank]->_traces.back().offset = int(v);
               Output_data[mpi_rank]->_traces.back().cdp = int(Input_it->cdp);    
               Output_data[mpi_rank]->_traces.back().scalco = short(fac1);
               it_left=Input_data->_traces.begin() +1;
           it_right=Input_data->_traces.end()-1;
    
    /// Loop over all samples: /////////////////////////////////////////////////
           for (int ii=0; ii<ns; ii++)
               {
               t0=ii*dtnew;
               Value=0.;
                   s_count=0;
            for (int is=0; is < width;is++){
                  sembl_up[is]=0.;sembl_down[is]=0.;
                }
            for ( Input_temp=it_left;Input_temp != it_right;Input_temp++) {
                            
                            y=0.5*(Input_temp->gy  + Input_temp->sy)*fac;                        /// ym difference
                            x=0.5*(Input_temp->gx  + Input_temp->sx)*fac;                        /// xm difference
                if ( (xm-x)*(xm-x)/(aperx*aperx) +  (ym-y)*(ym-y)/(apery*apery) > 1.) continue;        /// aperture
                tt=sqrt(t0*t0+4.*(xm-x)*(xm-x)/(v*v) +4.*(ym-y)*(ym-y)/(v*v));                         /// calculate tt hyperbola over all times 
                         idown=static_cast<int>(tt/dtnew);                                                      ///calculate the current index in arrays
                if (Input_temp > Input_data->_traces.end())  break;
                            ( (idown -rcoher ) > 0)? idown=idown-rcoher: idown=idown; 
                            ( (idown + rcoher) > ns)? idown=ns - rcoher: idown=idown;
                s_count++;                                                                /// trace counter
                for (int is=0; is<width;is++){                                                     /// semblance within the coherence band
                                Value_temp=float(Input_temp->data[idown+is]); 
                              sembl_up[is]+=Value_temp;
                              sembl_down[is]+=Value_temp*Value_temp;
                }
                            s=1.*s_count; 
                         }
             Val_up= Val_down=0.0;
                   for (int is=0; is<width;is++){
                 Val_up+=sembl_up[is]*sembl_up[is];
                  Val_down+=sembl_down[is];
               }
              Value=Val_up/Val_down/s; 
                if  ( isinf( Value) || isnan(Value ) )   Value=0.0;
                  Output_data[mpi_rank]->_traces.back().data[ii]=Value; 
             }
          }
       }
     }
    
    /// Write Data////////////////////////////////////////////////////////////////////////////////
    
    //    MPI_Barrier(MPI_COMM_WORLD);
    //    Output_data[mpi_rank]->write();
    
    int loopn;
      Output=new SUdata(Outputfile.c_str());
    
      for (int i=mpi_rank; i<mpi_size;i++) {
          Input_it =  Output_data[mpi_rank]->_traces.begin();
          loopn= Output_data[mpi_rank]->_traces.size();
    
          for (int i=1; i<loopn; i++) 
              {       
    
               Output->_traces.push_back( Trace(ns));
           Output->_traces.back().gx = int(Input_it->sx);
               Output->_traces.back().sx = int(Input_it->gx);
           Output->_traces.back().gy = int(Input_it->gy);
               Output->_traces.back().sy = int(Input_it->sy);
           Output->_traces.back().dt = int(Input_it->dt);
               Output->_traces.back().offset = int(Input_it->offset);
               Output->_traces.back().cdp = int(Input_it->cdp);    
               Output->_traces.back().scalco = short(Input_it->scalco);
             for (int ii=0; ii<ns; ii++) Output->_traces.back().data[ii]=Input_it->data[ii];
               Input_it++;
             }
       }
    
       MPI_Barrier(MPI_COMM_WORLD);
       Output->write();
    
    ///Time stop
     if (mpi_rank==0) { 
        endtime=MPI_Wtime();
        time=static_cast<int>(endtime-starttime);
        sec_to_time(time);
      }
    ///MPI end
      MPI_Barrier(MPI_COMM_WORLD);
      MPI_Finalize();
    }
    /// end main
    
    \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
    I shall be thank ful to you if anyone can help me. I am really facing troubles.

    smartkhawar@gmail.com
    Last edited by 2kaud; November 24th, 2016 at 06:46 AM. Reason: Code tags added

  2. #2
    2kaud's Avatar
    2kaud is offline Super Moderator Power Poster
    Join Date
    Dec 2012
    Location
    England
    Posts
    6,920

    Re: Azimuth search

    When posting code, please use code tags so that the code is readable. Go Advanced, select the formatted code and click '#'.

    I am working to make full azimuth search...I am really facing troubles
    So what is your specific question? What doesn't work, what specific 'troubles' are you having? What debugging of the code have you done?
    Last edited by 2kaud; November 24th, 2016 at 06:50 AM.
    All advice is offered in good faith only. All my code is tested (unless stated explicitly otherwise) with the latest version of Microsoft Visual Studio (using the supported features of the latest standard) and is offered as examples only - not as production quality. I cannot offer advice regarding any other c/c++ compiler/IDE or incompatibilities with VS. You are ultimately responsible for the effects of your programs and the integrity of the machines they run on. Anything I post, code snippets, advice, etc is licensed as Public Domain https://creativecommons.org/publicdomain/zero/1.0/ and can be used without reference or acknowledgement. Also note that I only provide advice and guidance via the forums - and not via private messages!

    C++17 Compiler: Microsoft VS2019 (16.4.0)

  3. #3
    Join Date
    Aug 2012
    Posts
    5

    Re: Azimuth search

    Dear Friend.
    Actually.
    I was only able to run this for single direction search. Trouble is that I am not so much expert in coding. I need the 180 degree azimuth search. This works. but only for the values it picks is in 0 degree. Please if you can help me out.

    waiting for help

    smartkhawar@gmail.com

  4. #4
    2kaud's Avatar
    2kaud is offline Super Moderator Power Poster
    Join Date
    Dec 2012
    Location
    England
    Posts
    6,920

    Re: Azimuth search

    I need the 180 degree azimuth search. This works. but only for the values it picks is in 0 degree.
    Its doubtful if anyone here is going to go through the code, understand how it works, fathom out how it's supposed to work and then change it.

    If part of the code is not working as expected, then use the debugger to trace through the code to see what the code is doing and to determine where it deviates from that expected from the design.

    Once you have a specific c++ question, then we'll be able to advise further.
    All advice is offered in good faith only. All my code is tested (unless stated explicitly otherwise) with the latest version of Microsoft Visual Studio (using the supported features of the latest standard) and is offered as examples only - not as production quality. I cannot offer advice regarding any other c/c++ compiler/IDE or incompatibilities with VS. You are ultimately responsible for the effects of your programs and the integrity of the machines they run on. Anything I post, code snippets, advice, etc is licensed as Public Domain https://creativecommons.org/publicdomain/zero/1.0/ and can be used without reference or acknowledgement. Also note that I only provide advice and guidance via the forums - and not via private messages!

    C++17 Compiler: Microsoft VS2019 (16.4.0)

  5. #5
    Join Date
    Aug 2012
    Posts
    5

    Re: Azimuth search

    Dear Sir Super Moderator., this issue is resolved by my self now. Please remove the thread. I will run the code now and if there is any C++ specific question i will ask you. Thank you . that you help people like me.
    Stay blessed. Greetings.

Posting Permissions

  • You may not post new threads
  • You may not post replies
  • You may not post attachments
  • You may not edit your posts
  •  


Windows Mobile Development Center


Click Here to Expand Forum to Full Width




On-Demand Webinars (sponsored)