123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417 |
- #include <iostream>
- #include <vector>
- #include <map>
- #include <utility>
- #include <numeric>
- #include <limits>
- #include <cmath>
- #include <TSystem.h>
- #include "filval/filval.hpp"
- #include "filval/root/filval.hpp"
- #include <boost/range/combine.hpp>
- #include <boost/format.hpp>
- #include "TrackingNtuple.h"
- #include "analysis/obj_types.cpp"
- using namespace std;
- using namespace fv;
- using namespace fv::root;
- typedef std::tuple<Seed, PixRecHit, Track, SimHit, SimTrack> MatchedTrack;
- typedef std::pair<std::vector<float>,std::vector<float>> pair_vec;
- #define HIT_TYPE_PIXEL 0
- #define HIT_TYPE_GLUED 1
- #define HIT_TYPE_STRIP 2
- #define PIXEL_BARREL 1
- #define PIXEL_ENDCAP 2
- std::map<int,string> subdet_names = {{1, "BPIX"}, {2, "FPIX"}};
- template<typename A, typename B>
- float displacement(const A& a, const B& b){
- return std::sqrt(pow(a.x-b.x, 2) +
- pow(a.y-b.y, 2) +
- pow(a.z-b.z, 2));
- }
- template<typename T>
- float rho(const T& t){
- return std::sqrt(pow(t.x,2)+pow(t.y,2));
- }
- bool in_det(const MatchedTrack &mt, int &&det, int &&layer){
- auto& hit = std::get<PixRecHit>(mt);
- return hit.det == det && hit.lay == layer;
- };
- float pseudorapidity(float z, float r){
- float theta = atan2(r, z);
- return -log(tan(theta/2.0));
- }
- float pseudorapidity(float x, float y, float z){
- float r = sqrt(x*x + y*y);
- return pseudorapidity(z, r);
- }
- template<typename T>
- float pseudorapidity(const T &t){
- return pseudorapidity(t.x, t.y, t.z);
- }
- vector<string> seedTypes =
- {"initialStepSeeds",
- "highPtTripletStepSeeds",
- "mixedTripletStepSeeds",
- "pixelLessStepSeeds",
- "tripletElectronSeeds",
- "pixelPairElectronSeeds",
- "stripPairElectronSeeds"};
- Value<vector<Seed>>* seeds;
- Value<vector<PixRecHit>>* pixrec_hits;
- Value<vector<Track>>* tracks;
- Value<vector<SimHit>>* sim_hits;
- Value<vector<SimTrack>>* sim_tracks;
- void register_objects(TrackingDataSet& tds){
- seeds = register_seeds(tds);
- pixrec_hits = register_pixrec_hits(tds);
- tracks = register_tracks(tds);
- sim_hits = register_sim_hits(tds);
- sim_tracks = register_sim_tracks(tds);
- }
- void setup_first_hit_pairs(TrackingDataSet& tds){
- typedef std::tuple<PixRecHit, SimHit> HitPair;
- auto& find_matched_nth_hit_in_layer =
- func<vector<HitPair>(vector<Seed>,
- vector<PixRecHit>,
- vector<SimHit>,
- vector<Track>,
- vector<SimTrack>,
- int, int, int)>("find_matched_nth_hit_in_layer",
- FUNC(([](const vector<Seed>& seeds,
- const vector<PixRecHit>& pixrec_hits,
- const vector<SimHit>& sim_hits,
- const vector<Track>& tracks,
- const vector<SimTrack>& sim_tracks,
- const int&& det,
- const int&& pix_layer,
- const int&& hit_number){
- vector<HitPair> matched_hits;
- for(const Track &trk : tracks){ // loop over all tracks
- const Seed &seed = seeds[trk.seedIdx];
- if(seed.hitIdx.size() <= hit_number) continue; // looking for hit_number'th hit, which this seed doesn't have
- if(seed.algoOriginal < 0 || seed.algoOriginal >= seedTypes.size()) continue;
- if(seed.hitType[hit_number] != HIT_TYPE_PIXEL) continue; // take only pixel hits for now
- const PixRecHit &rec_hit = pixrec_hits[seed.hitIdx[hit_number]];
- if(rec_hit.det == det && rec_hit.lay == pix_layer){
- // We have the RecHit we want to work with, now find a properly* matched SimHit
- if(rec_hit.simHitIdx.size() == 0) continue; // if rechit is matched to no simhits, give up.
- for(const int &simTrkIdx : trk.simTrkIdx){ // loop over SimTracks matched to Track
- for(const int& simHitIdx : sim_tracks[simTrkIdx].simHitIdx){ // loop over SimHits from SimTrack
- if(simHitIdx == rec_hit.simHitIdx[0]) // take first matched simhit (should be the closest one)
- matched_hits.push_back({rec_hit, sim_hits[rec_hit.simHitIdx[0]]});
- }
- }
- }
- }
- return matched_hits;
- })));
- auto barrel_val = constant("PIXEL_BARREL", PIXEL_BARREL);
- auto endcap_val = constant("PIXEL_ENDCAP", PIXEL_ENDCAP);
- auto first_hit = constant("1st", 0);
- auto second_hit = constant("2nd", 1);
- // First hits on inner three bpix layers
- auto first_hits_in_B1 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, barrel_val, constant("L1", 1), first_hit), "first_hits_in_B1");
- auto first_hits_in_B2 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, barrel_val, constant("L2", 2), first_hit), "first_hits_in_B2");
- // Second hits on outer three bpix layers
- auto second_hits_in_B2 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, barrel_val, constant("L2", 2), second_hit), "second_hits_in_B2");
- auto second_hits_in_B3 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, barrel_val, constant("L3", 3), second_hit), "second_hits_in_B3");
- auto second_hits_in_B4 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, barrel_val, constant("L4", 4), second_hit), "second_hits_in_B4");
- auto first_hits_in_F1 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, endcap_val, constant("L1", 1), first_hit), "first_hits_in_F1");
- auto first_hits_in_F2 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, endcap_val, constant("L2", 2), first_hit), "first_hits_in_F2");
- auto second_hits_in_F2 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, endcap_val, constant("L2", 2), second_hit), "second_hits_in_F2");
- auto second_hits_in_F3 = fv::apply(find_matched_nth_hit_in_layer,
- fv::tuple(seeds, pixrec_hits, sim_hits, tracks, sim_tracks, endcap_val, constant("L3", 3), second_hit), "second_hits_in_F3");
- // Even vs Odd Ladders
- auto even = constant("even", false);
- auto odd = constant("odd", true);
- auto& select_even_odd_ladder_blade_hit_pairs =
- func<vector<HitPair>(vector<HitPair>, bool)>("select_even_odd_ladder_blade_hit_pairs",
- FUNC(([](const vector<HitPair>& hit_pairs,
- const bool &&odd){
- vector<HitPair> even_pairs;
- for(const HitPair &hit_pair : hit_pairs){ // loop over all seeds
- if(std::get<PixRecHit>(hit_pair).ladder_blade % 2 == odd){
- even_pairs.push_back(hit_pair);
- }
- }
- return even_pairs;
- })));
- auto first_hits_in_B1_even_ladder = fv::apply(select_even_odd_ladder_blade_hit_pairs,
- fv::tuple(first_hits_in_B1, even), "first_hits_in_B1_even_ladder");
- auto first_hits_in_B1_odd_ladder = fv::apply(select_even_odd_ladder_blade_hit_pairs,
- fv::tuple(first_hits_in_B1, odd), "first_hits_in_B1_odd_ladder");
- auto first_hits_in_B2_even_ladder = fv::apply(select_even_odd_ladder_blade_hit_pairs,
- fv::tuple(first_hits_in_B2, even), "first_hits_in_B2_even_ladder");
- auto first_hits_in_B2_odd_ladder = fv::apply(select_even_odd_ladder_blade_hit_pairs,
- fv::tuple(first_hits_in_B2, odd), "first_hits_in_B2_odd_ladder");
- //Plots for dPhi of collections defined above
- auto& calc_dphi_v_eta = func<pair_vec(vector<HitPair>)>("calc_dphi_v_eta",
- FUNC(([](const vector<HitPair>& hit_pairs){
- vector<float> dphis;
- vector<float> etas;
- for(auto hit_pair : hit_pairs){
- auto& pixrec_hit = std::get<PixRecHit>(hit_pair);
- auto& sim_hit = std::get<SimHit>(hit_pair);
- dphis.push_back(atan2(sim_hit.x, sim_hit.y) - atan2(pixrec_hit.x, pixrec_hit.y));
- etas.push_back(pseudorapidity(pixrec_hit));
- }
- return std::make_pair(etas, dphis);
- })));
- TH2Params params_dphi = {"$\\eta$", 100, -4, 4,
- "$\\Delta \\phi$(rad)", 50, -.0015, .0015};
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B1",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B1)),
- "First Hit in BPIX-L1", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B2",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B2)),
- "First Hit in BPIX-L2", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B1_even_ladder",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B1_even_ladder)),
- "First Hit in BPIX-L1 - Even Ladders", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B1_odd_ladder",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B1_odd_ladder)),
- "First Hit in BPIX-L1 - Odd Ladders", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B2_even_ladder",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B2_even_ladder)),
- "First Hit in BPIX-L2 - Even Ladders", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_B2_odd_ladder",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_B2_odd_ladder)),
- "First Hit in BPIX-L2 - Odd Ladders", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_second_hits_in_B2",
- fv::apply(calc_dphi_v_eta, fv::tuple(second_hits_in_B2)),
- "Second Hit in BPIX-L2", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_second_hits_in_B3",
- fv::apply(calc_dphi_v_eta, fv::tuple(second_hits_in_B3)),
- "Second Hit in BPIX-L3", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_second_hits_in_B4",
- fv::apply(calc_dphi_v_eta, fv::tuple(second_hits_in_B4)),
- "Second Hit in BPIX-L4", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_F1",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_F1)),
- "First Hit in FPIX-L1", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_first_hits_in_F2",
- fv::apply(calc_dphi_v_eta, fv::tuple(first_hits_in_F2)),
- "First Hit in FPIX-L2", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_second_hits_in_F2",
- fv::apply(calc_dphi_v_eta, fv::tuple(second_hits_in_F2)),
- "Second Hit in FPIX-L2", params_dphi);
- tds.register_container<ContainerTH2Many<float>>("dphi_v_eta_second_hits_in_F3",
- fv::apply(calc_dphi_v_eta, fv::tuple(second_hits_in_F3)),
- "Second Hit in FPIX-L3", params_dphi);
- //Plots for dr of collections defined above
- auto& calc_dr_v_eta = func<pair_vec(vector<HitPair>)>("calc_dr_v_eta",
- FUNC(([](const vector<HitPair>& hit_pairs){
- vector<float> drs;
- vector<float> etas;
- for(auto hit_pair : hit_pairs){
- auto& pixrec_hit = std::get<PixRecHit>(hit_pair);
- auto& sim_hit = std::get<SimHit>(hit_pair);
- drs.push_back(displacement(pixrec_hit, sim_hit));
- etas.push_back(pseudorapidity(pixrec_hit));
- }
- return std::make_pair(etas, drs);
- })));
- TH2Params params_dr = {"$\\eta$", 100, -4, 4,
- "$\\Delta \\r$(rad)", 50, 0, .006};
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B1",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B1)),
- "First Hit in BPIX-L1", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B2",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B2)),
- "First Hit in BPIX-L2", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B1_even_ladder",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B1_even_ladder)),
- "First Hit in BPIX-L1 - Even Ladders", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B1_odd_ladder",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B1_odd_ladder)),
- "First Hit in BPIX-L1 - Odd Ladders", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B2_even_ladder",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B2_even_ladder)),
- "First Hit in BPIX-L2 - Even Ladders", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_B2_odd_ladder",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_B2_odd_ladder)),
- "First Hit in BPIX-L2 - Odd Ladders", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_second_hits_in_B2",
- fv::apply(calc_dr_v_eta, fv::tuple(second_hits_in_B2)),
- "Second Hit in BPIX-L2", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_second_hits_in_B3",
- fv::apply(calc_dr_v_eta, fv::tuple(second_hits_in_B3)),
- "Second Hit in BPIX-L3", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_second_hits_in_B4",
- fv::apply(calc_dr_v_eta, fv::tuple(second_hits_in_B4)),
- "Second Hit in BPIX-L4", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_F1",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_F1)),
- "First Hit in FPIX-L1", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_first_hits_in_F2",
- fv::apply(calc_dr_v_eta, fv::tuple(first_hits_in_F2)),
- "First Hit in FPIX-L2", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_second_hits_in_F2",
- fv::apply(calc_dr_v_eta, fv::tuple(second_hits_in_F2)),
- "Second Hit in FPIX-L2", params_dr);
- tds.register_container<ContainerTH2Many<float>>("dr_v_eta_second_hits_in_F3",
- fv::apply(calc_dr_v_eta, fv::tuple(second_hits_in_F3)),
- "Second Hit in FPIX-L3", params_dr);
- //Plots for dz of collections defined above
- auto& calc_dz_v_eta = func<pair_vec(vector<HitPair>)>("calc_dz_v_eta",
- FUNC(([](const vector<HitPair>& hit_pairs){
- vector<float> dzs;
- vector<float> etas;
- for(auto hit_pair : hit_pairs){
- auto& pixrec_hit = std::get<PixRecHit>(hit_pair);
- auto& sim_hit = std::get<SimHit>(hit_pair);
- dzs.push_back(sim_hit.z - pixrec_hit.z);
- etas.push_back(pseudorapidity(pixrec_hit));
- }
- return std::make_pair(etas, dzs);
- })));
- TH2Params params_dz = {"$\\eta$", 100, -4, 4,
- "$\\Delta z$(rad)", 50, -.01, .01};
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B1",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B1)),
- "First Hit in BPIX-L1", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B2",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B2)),
- "First Hit in BPIX-L2", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B1_even_ladder",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B1_even_ladder)),
- "First Hit in BPIX-L1 - Even Ladders", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B1_odd_ladder",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B1_odd_ladder)),
- "First Hit in BPIX-L1 - Odd Ladders", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B2_even_ladder",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B2_even_ladder)),
- "First Hit in BPIX-L2 - Even Ladders", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_first_hits_in_B2_odd_ladder",
- fv::apply(calc_dz_v_eta, fv::tuple(first_hits_in_B2_odd_ladder)),
- "First Hit in BPIX-L2 - Odd Ladders", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_second_hits_in_B2",
- fv::apply(calc_dz_v_eta, fv::tuple(second_hits_in_B2)),
- "Second Hit in BPIX-L2", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_second_hits_in_B3",
- fv::apply(calc_dz_v_eta, fv::tuple(second_hits_in_B3)),
- "Second Hit in BPIX-L3", params_dz);
- tds.register_container<ContainerTH2Many<float>>("dz_v_eta_second_hits_in_B4",
- fv::apply(calc_dz_v_eta, fv::tuple(second_hits_in_B4)),
- "Second Hit in BPIX-L4", params_dz);
- //Plots for drho of collections defined above
- auto& calc_drho_v_eta = func<pair_vec(vector<HitPair>)>("calc_drho_v_eta",
- FUNC(([](const vector<HitPair>& hit_pairs){
- vector<float> drhos;
- vector<float> etas;
- for(auto hit_pair : hit_pairs){
- auto& pixrec_hit = std::get<PixRecHit>(hit_pair);
- auto& sim_hit = std::get<SimHit>(hit_pair);
- drhos.push_back(rho(sim_hit) - rho(pixrec_hit));
- etas.push_back(pseudorapidity(pixrec_hit));
- }
- return std::make_pair(etas, drhos);
- })));
- TH2Params params_drho = {"$\\eta$", 100, -4, 4,
- "$\\Delta \\rho$(cm)", 50, -.01, .01};
- tds.register_container<ContainerTH2Many<float>>("drho_v_eta_first_hits_in_F1",
- fv::apply(calc_drho_v_eta, fv::tuple(first_hits_in_F1)),
- "First Hit in FPIX-L1", params_drho);
- tds.register_container<ContainerTH2Many<float>>("drho_v_eta_first_hits_in_F2",
- fv::apply(calc_drho_v_eta, fv::tuple(first_hits_in_F2)),
- "First Hit in FPIX-L2", params_drho);
- tds.register_container<ContainerTH2Many<float>>("drho_v_eta_second_hits_in_F2",
- fv::apply(calc_drho_v_eta, fv::tuple(second_hits_in_F2)),
- "Second Hit in FPIX-L2", params_drho);
- tds.register_container<ContainerTH2Many<float>>("drho_v_eta_second_hits_in_F3",
- fv::apply(calc_drho_v_eta, fv::tuple(second_hits_in_F3)),
- "Second Hit in FPIX-L3", params_drho);
- }
- void run_analysis(const vector<fv::util::DataFileDescriptor>& dfds, const string output_filename, bool silent){
- gSystem->Load("libfilval.so");
- auto replace_suffix = [](const std::string& input, const std::string& new_suffix){
- return input.substr(0, input.find_last_of(".")) + new_suffix;
- };
- string log_filename = replace_suffix(output_filename, ".log");
- fv::util::Log::init_logger(log_filename, fv::util::LogPriority::kLogInfo);
- TrackingDataSet tds(output_filename, dfds, "trackingNtuple/tree");
- /* tds.set_max_events(10); */
- register_objects(tds);
- setup_first_hit_pairs(tds);
- tds.process(silent);
- tds.save_all();
- }
- int main(int argc, char * argv[]){
- fv::util::ArgParser args(argc, argv);
- bool silent = args.cmdOptionExists("-s");
- string output_filename = args.cmdOptionExists("-o") ? args.getCmdOption("-o") : "output.root";
- if(args.cmdOptionExists("-F")){
- auto file_list = fv::util::read_input_list(args.getCmdOption("-F"));
- run_analysis(file_list, output_filename, silent);
- }else if(args.cmdOptionExists("-f")){
- string input_filename = args.getCmdOption("-f");
- string data_label = args.cmdOptionExists("-l") ? args.getCmdOption("-l") : "";
- fv::util::DataFileDescriptor dfd(input_filename, data_label);
- run_analysis(std::vector<fv::util::DataFileDescriptor>({dfd}), output_filename, silent);
- }else {
- cout << "Usage: ./tracking_validation (-s) {-o output_filename} -F datafiles.txt" << endl;
- cout << " ./tracking_validation (-s) {-l DATA_LABEL} {-o output_filename} -f treefile.root" << endl;
- }
- return 0;
- }
|