Changeset 2921 in Sophya

Timestamp:

Mar 27, 2006, 6:56:24 PM (20 years ago)

Author:

ansari

Message:

Modifications commandes FFT pour TD/TP SurvEnv - Reza 27/3/2006

File:

• trunk/SophyaPI/ProgPI/sopiamodule.cc (modified) (14 diffs)

trunk/SophyaPI/ProgPI/sopiamodule.cc

@@ -4 +4 @@
 #include <stdlib.h>
 #include <iostream>
+#include <ctype.h>
 #include <math.h>
 
@@ -32 +33 @@
 
 void SophyaFFT_forw(string& nom, string& nomout, string dopt);
-void SophyaFFT_back(string& nom, string& nomout, string dopt);
-void SophyaFFT_crefilter(string& nomvec, string& nomfilter, string func, string dopt);
-void SophyaFFT_filter(string& nom, string& nomfilter, string dopt);
+void SophyaFFT_back(string& nom, string& nomout, string dopt, bool forcez8=false);
+// void SophyaFFT_crefilter(string& nomvec, string& nomfilter, string func, string dopt);
+void SophyaFFT_filter(string& nom, string& filt, string& outvec, string dopt, bool fgvec=false);
 
 class sopiamoduleExecutor : public CmdExecutor {
@@ -196 +197 @@
 mpiac = omg.GetImgApp()->CmdInterpreter();
 
-string hgrp = "SophyaCmd";
 string kw,usage;
-
-kw = "powerspec";
-usage = "FFT on a vector -> Plots power spectrum ";
-usage += "\n Usage: powerspec vecName vecFFT [graphic_att] ";
-mpiac->RegisterCommand(kw, usage, this, hgrp);
+// ------------- Help group FFT 
+string hgrp = "FFT";
 
 kw = "fftforw";
-usage = "FFT on a vector -> Plots power spectrum (same as powerspec)";
-usage += "\n Usage: fftforw vecName vecFFT [graphic_att] ";
+usage = "FFT on a vector -> computes forward Fourier transform";
+usage += "\n Usage: fftforw vecSig vecSpec [graphic_att] ";
+usage += "\n vecSpec = FFTForward(vecSig) ";
+usage += "\n vecSig: Input data vector of type r_8 or complex<r_8> ";
+usage += "\n vecSpec: Output data vector of type complex<r_8> ";
+usage += "\n See also : fftback fftfilter fftfuncfilter ";
 mpiac->RegisterCommand(kw, usage, this, hgrp);
 
 kw = "fftback";
-usage = "FFT on a vector -> Plots power spectrum";
-usage += "\n Usage: fftback vecFFT vecName [graphic_att] ";
-usage += "\n  vecFFT is complex<r_8>";
-usage += "\n  If vecName already exist, it gives the type(r_8 or complex<r_8>)";
-usage += "\n  If vecName does not exist, return type is complex<r_8>";
-mpiac->RegisterCommand(kw, usage, this, hgrp);
-
-kw = "crefilter";
-usage = "Create a filter for vecFFT (vector complex<r_8>)";
-usage += "\n Usage: crefilter vecFFT filter f(x) [graphic_att] ";
+usage = "FFT on a vector -> computes backward Fourier transform ";
+usage += "\n Usage: fftback vecSpec vecSig [graphic_att] [C/Z] ";
+usage += "\n vecSig = FFTBackward(vecSpec) ";
+usage += "\n vecSpec: Input data vector of type complex<r_8> ";
+usage += "\n vecSig:  Output Data vector of type r_8 (default) ";
+usage += "\n    or complex<r_8> if C/Z specified or ";
+usage += "\n    vecSpec computed by fftforw on a complex vector";
+usage += "\n See also : fftforw fftfilter fftfuncfilter ";
 mpiac->RegisterCommand(kw, usage, this, hgrp);
 
 kw = "fftfilter";
-usage = "Filter (multiply) vecFFT (vector complex<r_8>) by filter (vector<r_8>)";
-usage += "\n Usage: fftfilter vecFFT filter [graphic_att] ";
-mpiac->RegisterCommand(kw, usage, this, hgrp);
+usage = "Filter (multiply) vecSpec (vector complex<r_8>) by Filter (vector<r_8>)";
+usage += "\n Usage: fftfilter vecSpec FilterVec vecFiltSpec [graphic_att] ";
+usage += "\n   vecFiltSpec(i) = vecSpec(i) * complex(FilterVec(i),0) ";
+usage += "\n See also : fftforw fftbackw fftfuncfilter ";
+mpiac->RegisterCommand(kw, usage, this, hgrp);
+
+kw = "fftfuncfilter";
+usage = "Filter (multiply) vecSpec (vector complex<r_8>) by FilterFunc(i) ";
+usage += "\n Usage: fftfilter vecSpec FilterFunc vecFiltSpec [graphic_att] ";
+usage += "\n   vecFiltSpec(i) = vecSpec(i) * complex(FilterFunc(i),0) ";
+usage += "\n See also : fftforw fftbackw fftfilter ";
+mpiac->RegisterCommand(kw, usage, this, hgrp);
+
+//--------------------------
+hgrp = "SophyaCmd";
 
 kw = "mollgridsph";
@@ -255 +266 @@
 if (kw == "powerspec" || kw == "fftforw") {
   if (tokens.size() < 2) {
-    cout << "Usage: powerspec nameVec vecFFT [graphic_att]" << endl;
+    cout << "Usage: fftforw vecSig vecSpec [graphic_att] " << endl;
     return(0);
     }
   if (tokens.size() < 3)  tokens.push_back((string)"");
   SophyaFFT_forw(tokens[0], tokens[1], tokens[2]);
-  }
+}
 else if (kw == "fftback") {
   if (tokens.size() < 2) {
-    cout << "Usage: fftback vecFFT nameVec [graphic_att]" << endl;
+    cout << "Usage: fftback vecSpec vecSig [graphic_att] [C/Z]" << endl;
     return(0);
     }
+  bool forcez8 = false;
+  if (tokens.size() > 3) {
+    char fcz8 = tolower(tokens[3][0]);
+    if ( (fcz8 == 'c') || (fcz8 == 'z') ) forcez8 = true;
+  }
   if (tokens.size() < 3)  tokens.push_back((string)"");
-  SophyaFFT_back(tokens[0], tokens[1], tokens[2]);
-  }
-else if (kw == "crefilter") {
+  SophyaFFT_back(tokens[0], tokens[1], tokens[2], forcez8);
+}
+else if ( (kw == "fftfilter") || (kw == "fftfuncfilter") ) {
   if (tokens.size() < 3) {
-    cout << "Usage: crefilter vecFFT filter f(x) [graphic_att]" << endl;
+    cout << "Usage: fftfilter/fftfunfilter vecSpec Filter vecFiltSpec [graphic_att]" << endl;
     return(0);
     }
   if (tokens.size() < 4)  tokens.push_back((string)"");
-  SophyaFFT_crefilter(tokens[0], tokens[1], tokens[2], tokens[3]);
-  }
-else if (kw == "fftfilter") {
-  if (tokens.size() < 2) {
-    cout << "Usage: fftfilter vecFFT filter [graphic_att]" << endl;
-    return(0);
-    }
-  if (tokens.size() < 3)  tokens.push_back((string)"");
-  SophyaFFT_filter(tokens[0], tokens[1], tokens[2]);
-  }
+  bool fgvec = false;
+  if (kw == "fftfilter") fgvec = true; 
+  SophyaFFT_filter(tokens[0], tokens[1], tokens[2], tokens[3], fgvec);
+}
 
 else if (kw == "prjmoldefval") {
@@ -292 +302 @@
     }
   Set_Project_Mol_DefVal(atof(tokens[0].c_str()) );
-  }
+}
+
 else if ( (kw == "mollgridsph") || (kw == "mollgrid") ) {
    NamedObjMgr omg;
@@ -394 +405 @@
   ffts.FFTForward(*vin, *vout);
   cout<<"...Output TVector<complex<r_8>> of size: "<<vout->NElts()<<endl;
+  vout->Info()["FFTTYPE"]="F:R8->Z8";
 } else {
   cout<<"SophyaFFT_forw() - Computing FFT of vector TVector<complex<r_8>> of size "<<vinc->NElts()<< endl;
   ffts.FFTForward(*vinc, *vout);
   cout<<"...Output TVector<complex<r_8>> of size: "<<vout->NElts()<<endl;
+  vout->Info()["FFTTYPE"] = "F:Z8->Z8";
 }
 // cout << " SophyaFFT - FFT done " << endl;
@@ -410 +423 @@
 
 /*  Nouvelle-Fonction */
-void SophyaFFT_back(string& nom, string& nomout, string dopt)
+void SophyaFFT_back(string& nom, string& nomout, string dopt, bool forcez8)
 {
   //Timer tm("powerspec");
@@ -428 +441 @@
 TVector< r_8 > * vout = NULL;
 TVector< complex<r_8> > * voutc = NULL;
-obj=omg.GetObj(nomout);
-if(obj == NULL) {
+
+if (forcez8 || (vinc->Info().GetS("FFTTYPE") == "F:Z8->Z8") ) {
   voutc = new TVector< complex<r_8> >;
-} else {
-  vout = dynamic_cast<TVector<r_8>*>(obj);
-  if(vout==NULL) voutc = dynamic_cast<TVector< complex<r_8> >*>(obj);
-}
-if(vout==NULL && voutc==NULL) {
- cout << "SophyaFFT_back() Error , Output objet n'est pas un vecteur r_8 ou complex<r_8>" <<  endl;
- return;
+} 
+else {
+  vout = new TVector< r_8 >;
 }
 
@@ -445 +454 @@
   ffts.FFTBackward(*vinc, *vout);
   cout<<"...Output TVector<r_8> of size: "<<vout->NElts()<<endl;
+  vout->Info()["FFTTYPE"]="B:Z8->R8";
   omg.AddObj(vout, nomout);
 } else {
   ffts.FFTBackward(*vinc, *voutc);
   cout<<"...Output TVector<complex<r_8>> of size: "<<voutc->NElts()<<endl;
+  voutc->Info()["FFTTYPE"]="B:Z8->Z8";
   omg.AddObj(voutc, nomout);
 }
@@ -471 +482 @@
   return;
 }
-int n = vfft->Size();
+sa_size_t n = vfft->Size();
 if(n<=0) return;
 
@@ -479 +490 @@
 
 /*  Nouvelle-Fonction */
-void SophyaFFT_filter(string& nom, string& nomfilter, string dopt)
+void SophyaFFT_filter(string& nom, string& filt, string& outvec, string dopt, bool fgvec)
 {
 NamedObjMgr omg;
@@ -487 +498 @@
   return;
 }
-
+TVector< complex<r_8> >* vfft = dynamic_cast<TVector< complex<r_8> >*>(obj);
+if(vfft == NULL) {
+  cout<<"SophyaFFT_filter() Error , Objet "<<nom<<" n'est pas un Tvector<complex<r_8>>"<< endl;
+  return;
+}
+int n = vfft->Size();
+if(n<=0) {
+  cout<<"SophyaFFT_filter() Error , vector size = 0 " << endl;
+  return;
+}
+
+string nomfilter = "/autoc/fft_filter";
+if (fgvec) {
+  nomfilter = filt; 
+  cout<<"SophyaFFT_filter(): using vector " << filt << " for filtering spectrum ..." << endl;   
+}
+else {
+  cout<<"SophyaFFT_filter(): Creating filter vector "<<nomfilter
+      <<" from function " << filt << endl; 
+  string fcrdopt = "nodisp";
+  omg.GetServiceObj()->PlotFunc(filt, nomfilter, 0., (double)n, n, fcrdopt); 
+}
 AnyDataObj* objfilter=omg.GetObj(nomfilter);
 if(objfilter == NULL) {
@@ -494 +526 @@
 }
 
-TVector< complex<r_8> >* vfft = dynamic_cast<TVector< complex<r_8> >*>(obj);
-if(vfft == NULL) {
-  cout<<"SophyaFFT_filter() Error , Objet "<<nom<<" n'est pas un Tvector<complex<r_8>>"<< endl;
-  return;
-}
-
 TVector<r_8>* vfilter = dynamic_cast<TVector<r_8>*>(objfilter);
 if(vfft == NULL) {
@@ -506 +532 @@
 }
 
-int n = (vfilter->Size()<vfft->Size()) ? vfilter->Size(): vfft->Size();
-cout<<"...filtering from 0 to "<<n<<endl;
-for(int i=0;i<n;i++) (*vfft)(i) *= (*vfilter)(i);
-
-if(dopt.size()>0) omg.DisplayObj(nom,dopt);
+n = (vfilter->Size()<vfft->Size()) ? vfilter->Size(): vfft->Size();
+
+cout<<" SophyaFFT_filter(): filtering from 0 to "<<n<<endl;
+TVector< complex<r_8> > * filtfft = new TVector< complex<r_8> >(*vfft, false);
+for(int i=0;i<n;i++) (*filtfft)(i) *= (*vfilter)(i);
+
+cout<<"SophyaFFT_filter(): Adding obj: " << outvec << endl;   
+omg.AddObj(filtfft, outvec);
+omg.DisplayObj(outvec, dopt);
 return;
 }