packages.ufo: add OutputRawDataFile

packages/ufo/src/plot.cpp

@@ -307,8 +307,10 @@ void ufo::plot_point(std::string config_file)
     std::optional<std::vector<std::pair<double, std::string>>> XTicks;
     // 是否输出图片
     std::optional<std::string> OutputPictureFile;
-    // 是否输出数据，可以进一步使用 matplotlib 画图
+    // 是否输出插值后数据，可以进一步使用 matplotlib 画图
     std::optional<std::string> OutputDataFile;
+    // 是否输出插值前数据，可以配合 phonopy 结果深入研究
+    std::optional<std::string> OutputRawDataFile;
   };
 
   // 根据 q 点路径, 搜索要使用的 q 点，返回的是 q 点在 QpointData 中的索引
@@ -346,8 +348,9 @@ void ufo::plot_point(std::string config_file)
     std::vector<double> result(resolution);
     for (auto& mode : qpoint.ModeData)
     {
-      std::ptrdiff_t index = (mode.Frequency - frequency_range[0]) / (frequency_range[1] - frequency_range[0])
-        * resolution;
+      double index_double = (mode.Frequency - frequency_range[0]) / (frequency_range[1] - frequency_range[0])
+        * (resolution - 1);
+      std::ptrdiff_t index = std::round(index_double);
       if (index >= 0 && index < static_cast<std::ptrdiff_t>(resolution)) result[index] += mode.Weight;
     }
     return log.rtn(result);
@@ -364,7 +367,7 @@ void ufo::plot_point(std::string config_file)
     biu::Logger::Guard log;
     auto f = matplot::figure(true);
     auto ax = f->current_axes();
-    auto image = ax->plot(values);
+    auto image = ax->area(values, 0, false, "");
     ax->y_axis().reverse(false);
     ax->x_axis().tick_values(x_ticks);
     ax->x_axis().tick_length(1);
@@ -423,4 +426,6 @@ void ufo::plot_point(std::string config_file)
       .write("XTickLabels", x_ticklabels)
       .write("InterpolationResolution", input.InterpolationResolution)
       .write("FrequencyRange", input.FrequencyRange);
+  if (input.OutputRawDataFile)
+    std::ofstream(*input.OutputRawDataFile) << YAML::Node(unfolded_data.QpointData[qpoint_index]);
 }